Atul D Sontakke
email: atul.sontakke@chimie-paristech.fr

Jean-Marie Mouesca

Victor Castaing

Alban Ferrier

Mathieu Salau

Isabelle Gautier-Luneau

Vincent Maurel

Alain Ibanez

Bruno Viana
email: bruno.viana@chimie-paristech.fr

Time-gated triplet-state optical spectroscopy to decipher organic luminophores embedded in rigid matrices †

Wet-chemically synthesized inorganic materials often exhibit luminescence behavior. We have recently shown that the amorphous yttrium-aluminium-borate (a-YAB) powders obtained by sol-gel and modified Pechini methods exhibit organic impurities, responsible for their intense visible photoluminescence and phosphorescence afterglow. However, the heterogeneity of impurity organic compounds and difficulties in their intact extraction from the solid inorganic host matrix limit the extraction-based chemical analysis for luminophore identification. Here, we propose a photo-physical route based on time-gated tripletstate optical spectroscopy (TGTSS) to construct the electronic structures of the trapped unknown luminophores, which successfully illustrates the luminescence properties of a-YAB powders in more detail and also provides important insights intrinsic to the nature of the luminophores. The experimental results accompanied with TD-DFT calculations of the theoretical electronic structures thus help us to propose the probable luminophore compounds trapped in rigid a-YAB matrices. We anticipate that the present approach will open new opportunities for analyzing similar complex luminescent materials, including carbon dots, graphene oxides, etc., which is vital for their improvement.

Introduction

Organic compounds in solution and vapor phases are rather easy to distinguish using state-of-the-art analytical techniques even in trace amounts. In rigid matrices, however, the molecular heterogeneity added with its complex bonding within rigid organic-inorganic frameworks unfavorably reduces the detection sensitivity, which then depends on their extraction yield from the host matrix. 1,2 Chemical extraction is highly controlled by the matrix properties, such as its ability to get dissolved in organic solvents or corrosive reactants (mild acids, bases, etc.). In the case of inorganic matrices, the requirements of relatively harsh reactants complicate the extraction of organic compounds in their original form. 3 This reduces the reliability of extraction based analytical techniques, thereby restricting their analysis to less precise solid-state methods. The analyses become even more challenging if one has to deal with particular functional molecules amongst heterogeneous molecular crowd strongly confined in rigid host framework.

We experienced similar difficulty in pinpointing the luminescent centers trapped in amorphous Yttrium-Aluminium-Borate (a-YAB) matrices. These a-YAB micrometer-sized powders prepared by modified Pechini and sol-gel methods give intense visible emissions with potential for developing metal ion-free luminescent materials, particularly phosphors for solid-state lighting. [4][5][6] Such intense visible luminescence has also been reported in many other inorganic materials obtained by various wet-chemical routes. [7][8][9][10][11] During the synthesis process, the metal-organic precursors undergo a complex chain of reactions, including polyesterification (Pechini route) or polycondensation (sol-gel chemistry) reactions for the preparation of initial resins or gels, followed by degradation and the associated structural transformation during the thermal treatments performed under controlled atmospheric conditions (pyrolysis or calcination). This leads to a variety of organic clusters trapped inside the solid inorganic matrix. Indeed, we showed by differential thermal analysis and thermogravimetry, coupled with mass spectrometry, that the a-YAB powders contain significant amounts of carbonaceous (organic) compounds (200-10 000 ppm) arising from precursors and solvents. 4,5 Moreover, we were able to make progress in establishing a correlation of photoluminescence (PL) emission with the presence of aromatic compounds such as small polycyclic aromatic hydrocarbon (PAH) molecules based on the observations of solid-state NMR, EPR-ENDOR and associated DFT calculations. 4 However, still precise information on luminescent compounds is not fully at reach, which is crucial in improving the performance of these promising metal ion-free luminescent materials.

Besides analytical methods, the photo-physical studies can be of great use to resolve the electronic structure and luminescence mechanism of the luminophores. The steady-state PL reveals excitation-emission characteristics, which is informative to understand the preliminary electronic structure of the luminophores. Unfortunately, a vast amount of organic compounds are luminescent in nature and differentiating them based on mere PL-PLE (PL excitation) or decay characteristics is simply impossible. Nevertheless, additional details about the location of singlet (S i ) and triplet (T i ) states, the localized or delocalized nature of excited states, the singlet-triplet energy difference (DE ST ), excited state dynamics, etc., could certainly narrow down the search and enrich the understanding of the trapped luminescent centers. This additional information provides better prospects for experimental and/or theoretical studies tallying with known probable luminophore candidates. Time-gated spectroscopy reveals the triplet-state dynamics, 12 which coupled with steady-state spectroscopy, can allow constructing detailed electronic structure and obtaining important insights on the excited state dynamics of the luminophores. Moreover, time-gated spectroscopy unveils weak effects, which are otherwise masked due to strong fluorescence in steady-state measurements. The triplet-state efficiency is generally weak at ambient temperature, but in rigid matrices, it is rather stabilized and yields efficient phosphorescence even at room temperature. 13 In a-YAB powders, the afterglow emission lasts for several seconds at room temperature (RT). 6

Results and discussion

Fig. 1a shows the steady-state PL-PLE contour plot of SG-450 a-YAB powder, synthesized from a sol-gel resin calcined at 450 1C. The powder exhibits blue PL emission corresponding to two distinct excitation contours peaking at about 255 nm and 350 nm. A careful observation reveals that the 350 nm PLE contour is non-symmetric and the peak maximum points to a relatively shorter PL wavelength than for the 255 nm PLE contour. This represents an important revelation about the emitting centers and their mechanism, however it is not possible to resolve it more from the steady-state luminescence properties. To gain better insights, the time-gated afterglow emission (PhL) was recorded and presented as the PhL-PhLE (PhL excitation) contour plot in Fig. 1b. Unlike the steady-state, the afterglow emission contour exhibits intense 255 nm excitation monitoring the blue emission but the corresponding 350 nm excitation is missing. This suggests that only high-energy excitation is contributing to the afterglow emission. Moreover, additional excitation contours are present corresponding to afterglow emissions at longer wavelengths with decreasing strengths. These weak excitation contours peaking at about 295 nm, 325 nm, 345 nm, etc., suggest that the a-YAB powder exhibits multiple luminescent centers coexisting at different proportions and are denoted as a, b, g, d, etc., as in Fig. 1b.

Among the different luminescent centers, the a ones are most dominant in SG-450 powder, followed by b, g, d, etc., with decreasing weights. Since the additional luminophores' transitions overlap the lower energy excitation of a centers, the PLE contour exhibits a non-symmetrical nature in this region (Fig. 1a). Moreover, this also shifts the contour peak position, which therefore represents the convoluted peak due to multiple contributions. Monitoring the PL higher energy edge region at 380-400 nm in Fig. 1a, which is less influenced by additional longer wavelength emitting luminophores, reveals the actual PLE positions of a emitting centers. Accordingly, the a centers exhibit excitations at 255 nm (S 2 ) and 330 nm (S 1 ).

Another interesting observation that can be made from Fig. 1 is that the steady-state PL and afterglow emission of SG-450 powder under 255 nm excitation exhibit similar profiles but differ when excited at lower energy excitations (Fig. 2a). 6 It is critical to correctly analyze the PL and afterglow emission in luminescent centers, which otherwise mislead the origin of the afterglow mechanism. In organic luminophores, the slow afterglow can be caused by the spin-forbidden triplet-to-singlet (T 1 -S 0 ) transition, thermally-assisted delayed fluorescence (TADF) or inverted singlet-triplet (iST), also known as triplettriplet annihilation mechanisms (ESI, † Scheme S1). [13][14][15] Among them, the first mechanism represents a pure phosphorescence transition from the triplet excited state to the singlet ground state, whereas the latter two mechanisms are effectively singletsinglet (S 1 -S 0 ) transitions delayed by reverse intersystem crossing (rISC). This incorporates a fundamental difference to distinguish these processes experimentally. The phosphorescence (T 1 -S 0 ) transition generally occurs at a longer wavelength with respect to the fluorescence (S 1 -S 0 ) transition following the classical Hund's rule, but the TADF and iST mechanisms exhibit the same fluorescence and delayed afterglow emission profiles owing to their identical origin (S 1 -S 0 ). 16 Interestingly, in SG-450, the PL and afterglow emission profiles look similar under 255 nm excitation (Fig. 2a), thus hinting the presence of the TADF or iST process. However, it is to be noted that the steady-state PL is composed of both prompt fluorescence and slow afterglow emission. The latter is significantly high at RT under 255 nm excitation, while it is rather weak under 350 nm excitation (Fig. 2b).

To reveal the peak position of the fluorescence transition, it is necessary to get rid of the afterglow emission in PL signals. This can be achieved in the PL spectra recorded at a higher temperature, as it assists in thermally accelerated diffusion from the excited triplet state to quenching sites. 17 Fig. 2c shows the temperature-dependent PL spectra of SG-450 powder. At 12 K, the PL peaks at 435 nm but gradually shifts towards a shorter wavelength with the increase in temperature. This is due to the fact that the PL at low temperature is dominated with slow afterglow emission, but its contribution rapidly vanishes with the increase of temperature, thereby the fluorescence peak resolves more prominently. This is also evident from the PL thermal stability plot (Fig. 2d), exhibiting step-like thermal quenching behavior. The initial continuous decrease in PL up to 400 K is due to the thermal quenching of afterglow emission contribution. 16 Above 400 K, the PL attains a relatively stable plateau representing the thermal stability of fluorescence emission. At 600 K, the PL exhibits a prominent peak at 410 nm, suggesting the peak wavelength of the fluorescence transition of a luminescent centers. Considering the thermally influenced red shift, the fluorescence peak at 300 K may be located at a slightly lower wavelength. This confirms that the slow afterglow is significantly red shifted with respect to fluorescence and is due to the phosphorescence mechanism. The results suggest that the fluorescence (B405 nm) and phosphorescence (B440 nm) transitions in the a centers are separated by about 0.25 eV, which may account for the singlet-triplet energy difference (DE ST ). Ideally, DE ST represents the S 1 -T 1 energy splitting observed through the zero phonon line (0-0) transitions from the singlet (S 1 -S 0 ) and triplet (T 1 -S 0 ) states. However, this is rather difficult to measure accurately due to their different degrees of symmetry forbiddenness. 18 Moreover, in the amorphous solid host, the random variation of dielectric polarization of the surrounding neighborhood leads to a Gaussian-type density of state distribution even at low temperature, unlike in the gas or liquid phase where the vibronic sidebands are vividly distinguishable. 17,19 Therefore, the peak positions of the fluorescence and phosphorescence transitions are used to estimate the DE ST value in the present case.

Though the triplet state emission is at a lower energy than the singlet state emission, the triplet state gets populated only by 255 nm excitation, whereas both 255 nm and 330 nm excitation bands contribute to the singlet emission in the a centers. This hints that the intersystem crossing (ISC) coupling interactions of singlet excitation states with the triplet state may not be the same. Fig. 3 shows the diffuse reflectance spectrum of SG-450 powder along with its steady-state PL excitation (PLE) and phosphorescence excitation (PhLE) monitoring the emission peak. The diffuse reflectance spectrum is plotted in the absorption mode for comparison. This shows more resemblance with the phosphorescence excitation spectrum, where the absorption maximum is at about 250 nm and gradually decreases for longer wavelengths. As the absorption represents the relative strength of transitions, it can be concluded that the oscillator strength of S 2 transitions contributing to both PL and phosphorescence is stronger than the S 1 PLE transition contributing only to the steady-state PL. The high oscillator strength coupled with ease of populating the forbidden triplet state indicates that the S 2 state may be delocalized in nature and exhibits a lower binding energy for excited electrons to allow spin flips, thereby favoring an efficient intersystem crossing (ISC) probability to the triplet state. Such properties are intrinsic to charge transfer (CT) or excitonic transitions. [13][14][15][16]20 This suggests that the 255 nm (S 2 ) excitation of the a centers is of delocalized CT nature. Similarly, the 295 nm, 325 nm, 345 nm, etc. excitation bands contributing to phosphorescence in b, g, d, etc., centers are delocalized S 2 states, respectively. This further illustrates why the phosphorescence in SG-450 is significantly red shifted with respect to PL under 350 nm excitation and possesses a relatively weak afterglow (Fig. 1 and2a, b). Accordingly under 350 nm excitation, the steady-state PL originates from dominant a centers accompanied by small contributions from additional overlapping centers, but the phosphorescence mainly arises from the d ones (with small fractions from g and e).

The above results are compiled to construct an energy-level diagram representing the electronic structure of a luminescent centers (Fig. 4). The other luminophores (b, g, d, etc.) do possess a similar energy level scheme but with lower energy transitions. It is evident that the S 2 state, upon excitation, can populate both S 1 and T 1 states by internal conversion (IC) and ISC mechanisms, giving rise to fluorescence and phosphorescence, respectively. In contrast, the ISC from the S 1 state to T 1 is very weak, thereby only contributing to fluorescence. Moreover, the T 1 to S 1 rISC is also missing, even DE ST is only about 0.25 eV. Normally, DE ST o 0.3 eV leads to the TADF process, but in the present case, no TADF is observed in the afterglow emission. This may be due to the fact that the TADF process also requires free dynamical movement between donor (D)-acceptor (A) molecular groups, 14,21 which might be constrained in the a-YAB matrix owing to the rigid inorganic glassy encapsulation. This suggests that the restricted rISC probability accompanied by low nonradiative losses due to the strong confinement of luminescent centers in rigid matrices is primarily responsible to increase the triplet state excitation lifetime that yields ultra-long phosphorescence even at elevated temperatures.

In an earlier report, we identified the presence of certain organic compounds confined in the a-YAB matrix using solidstate NMR and EPR spectroscopies. 4 The pulsed-ENDOR associated DFT calculations revealed that the small (2-3 rings) PAH molecules exhibit good correlation with the PL properties of a-YAB powders. 4 The electronic structures and S-S/S-T excitation energies of these small PAH molecules are calculated here using TD-DFT methods (ESI, † Table S1). Among the several selected compounds, the electronic structure of the phenalene molecule possesses a significant resemblance with the allowed transition of a luminophores presented in Fig. 4. However, the first theoretical triplet-state energy of phenalene is lower than the experimental results. It is to be noted that the calculated molecules are just the simplest working models (based on the radical species identified by magnetic resonance techniques), 4 while the true trapped molecule(s) will most probably be chemically related, exhibiting additional pending functional group(s), which may cause the favorable lift in the triplet-state energy. 21 The solid host matrix due to its polarization influence also brings a marked reduction in DE ST , associated with the hybrid charge transfer type characters of the singlet-triplet energy levels. 22,23 The important point here is that phenalene (or chemically-related molecules) satisfies both magnetic (NMR/EPR) resonance measurements and the absorption and PL-PLE experiments.

The evolution of aromatic NMR and ENDOR signals particularly followed the broadening of PL profiles in a-YAB powders synthesized at different temperatures. 4 Fig. 5a shows the steady-state PL of a-YAB powders calcined in a temperature range of 450-600 1C. At 450 1C, SG-450 powder exhibits intense blue emission, but the emission becomes broader for higher calcination temperatures and shifts towards the white color with gradual appearance followed by strengthening of longer wavelength emissions. This behavior is related to the evolution and relative strengths of the emitting centers with calcination temperature. In SG-450, the a centers are most prominent leading to intense violet-blue emission. For the powders calcined at higher temperatures, the dominance systematically shifts towards longer wavelength emitting centers, thereby incorporating a gradual change in emission spectrum and PL color coordinates (Fig. 5b). This becomes clear from the PL-PLE contour plots in Fig. 5c andd. The SG-470 sample exhibits distinct contours similar to SG-450 (Fig. 1a), but with more asymmetry in the lower energy region contour. This is due to the fact that the b and other additional centers (g, d, e, etc.) are relatively more prominent in SG-470 over SG-450 powder, although the a luminophores are still the most dominant ones. This systematic admixture of longer wavelength emitting centers with the increase of calcination temperature up to 600 1C induces the broadening in PL spectra (Fig. 5d; and ESI, † Fig. S1 andS2).

The relative strengthening of the longer wavelength emitting center could be confirmed from the time-gated phosphorescence plot of higher calcination temperature powders. Fig. 6a shows the phosphorescence contour plot of SG-500, still possessing the dominant contribution from a centers' S 2 state, but the relative contribution from b, g, d, e, etc., also increased notably compared to SG-450 (Fig. 1b). This systematic increase in longer wavelength emitting centers' fraction is more vivid in SG-570 (Fig. 6b). It can also be noticed that the phosphorescence profile is broadened and red shifted in SG-500 and SG-570 powders over SG-450 powder (Fig. 6c). This is due to the radiation trapping effect as the longer wavelength emitting centers become more abundant with higher calcination temperature and their absorption overlaps with the emission of shorter wavelength emitting centers as seen from the allowed energy level schemes of the first few prominently deducible emitting centers (a, b, g, d and e) in a-YAB powders (Fig. 7). The transition energies were obtained from the PL-PLE measurements of various samples. The relative prominence of longer wavelength emitting centers is also reflected in the absorption profiles of the a-YAB powders, which exhibit a gradual broadening as a function of calcination temperature (ESI, † Fig. S3).

The systematic evolution in PL allows smooth tuning of PL colors from blue to yellow intersecting with white color chromaticity. This is, in particular, an interesting property in view of phosphors for solid-state lighting devices with tunable chromaticity. Moreover, the continuous distribution of PL emissions covering the whole visible range helps to achieve high color rendering ability and reasonable luminous efficacy. 24 The color rendering index, R a , the correlated color temperature, CCT, the luminous efficacy of radiation (LER) and the luminescence internal quantum efficiency (iQY) values are presented in Table 1. Both R a and LER are maximum for the powders calcined at 530-570 1C, which are also optimized powders in view of relative intensity (inset Fig. 5a) and white color chromaticity (CCT B 5500 K; Fig. 5b). The iQY is about 30% in the powders used, which suggests that a careful tailoring of small PAH molecules trapped in aluminium borate matrices can potentially lead to efficient phosphors.

Conclusions

In summary, time-gated triplet-state spectroscopy has been used to construct the electronic structures of organic luminophores trapped in the a-YAB glassy matrix (micrometer-sized powders). The results successfully explained the luminescence mechanism in a-YAB powders in more detail. TD-DFT analysis revealed that among directly identified PAH compounds via magnetic resonance techniques, the phenalene-related compounds exhibit most resemblance with the experimentally derived electronic structure of the luminophores in the a-YAB matrix. The results demonstrated that time-gated spectroscopy is a vital technique to understand the complex PL mechanism in such defect-induced phosphors and to specify the nature of the confined molecular impurities in the rigid network, which play the role of PL emitting centers. method. In the first solution, 2.7 mL of B(OC 2 H 5 ) 3 (99%, Sigma Aldrich) was added in 3.5 mL of absolute ethanol. The addition of 3.1 mL Al(OCH(CH 3 )CH 2 CH 3 ) 3 (98%, STREM Chemicals) at 50 1C gives rise to the transesterification reaction among the reactants resulting in the formation of the aluminium borate network. In the second solution, 1.532 g Y(NO 3 ) 3 Á6H 2 O (99.9%, STREM Chemicals) was dissolved in 0.9 mL of propionic acid (CH 3 CH 2 COOH, 99%, Sigma Aldrich) at 50 1C. This solution was then added to the first one, inducing the dissolution of the aluminium borate solid and leading to a translucent yellowish solution. The molar ratio of Y : Al : B : propionic acid is equal to 1 : 3 : 4 : 3. The resulting solution was then placed under reflux (2 h) to maximize the hydrolysis and condensation reactions of the reagents. Then, it was evaporated at 70 1C to obtain a dry resin. This resin was carefully hand-milled before thermal treatments at various calcination temperatures, between 450 and 600 1C, under an O 2 atmosphere. The heating rate was kept at 2 1C min À1 and the soaking duration was one hour.

Experimental

(b) Experimental methods

The PL measurement was carried out using a broadband Xe plasma lamp (Energetiq EQ-99X) attached with a Jobin-Yvon UV monochromator for excitation, and the signals were recorded using an Avantes CCD spectrophotometer (AvaSpec-2048 TEC). The spectrometer was calibrated for a relative spectral power distribution using a halogen lamp (Ava-Hal-Cal) and the excitation lamp intensity correction was made using a UV-enhanced Si photo power meter (Thorlabs PM100D). The same setup was used for PL-PLE contour mapping, where the PL was recorded for a series of excitation wavelengths with 5 nm step resolution.

In the case of afterglow contour mapping (PhL-PhLE), the setup including the broadband Xe plasma lamp attached with the monochromator was used for excitation and the signals were recorded using an Acton SpectraPro 2150 Dual Grating Monochromator and a Pixis 100 CCD camera with a manual shutter to cease the excitation. The PhL-PhLE contour represents time-gated phosphorescence signals at 1 s (0.5 s in the case of SG-570) after ceasing the excitation. For the afterglow decay, the same setup was used as above for excitation and emission. The PL thermal stability measurements were carried out using a Sumitomo Cryogenics HC-4E closed cycle He cryostat connected with a Lakeshore 340 Temperature controller. The PL was recorded under continuous excitation of a Xe lamp during controlled heating from 10 K to 600 K. The UV-visible-NIR reflectance spectra of the powder samples were recorded using a Cary-Eclipse UV-visible-NIR spectrophotometer in diffuse reflectance mode.

(c) TD-DFT analysis

All DFT calculations have been performed with the ADF (Amsterdam Density Functional) code developed by E. J. Baerends and co-workers 25 using triple-zeta (TZ2P) basis sets (no frozen core). The S-S electronic transitions, oscillator strengths, and S-T transitions are computed using the B3LYP 26 exchangecorrelation potential (20% Hartree-Fock) and are listed in the ESI † (Table S1).

Fig. 1

 1 Fig. 1 Steady-state (a) and time-gated (b) luminescence contour maps of SG-450 a-YAB powder. In (a), the diagonal features are excitation pump signals. In the time-gated plot (b), data represent the emission signals at 1 s after ceasing the excitation. The intensity scale bar is presented on the right side.

Fig. 3

 3 Fig. 3 Absorption spectrum of SG-450 powder (a); the excitation spectrum for phosphorescence (PhLE) (b); and steady-state (PLE) emission (c) monitoring 440 nm emission.

Fig. 2

 2 Fig. 2 Steady-state (PL) and time-gated afterglow (PhL) emission spectral profiles of SG-450 powder under 255 nm and 350 nm excitations (a); the corresponding afterglow decay kinetics normalized with respect to the steady-state PL (b); temperature dependence of the steady-state PL spectra (c); and the variation of the integrated PL intensity as a function of temperature (d).

Fig. 4

 4 Fig. 4 Jablonski diagram representing the electronic transitions of a centers in SG-450 powder.

Fig. 5

 5 Fig. 5 Steady-state PL spectra of a-YAB powders calcined in a temperature range of 450-600 1C recorded under 365 nm excitation (a). The inset indicates the integrated PL intensity (radiance) and the luminance contribution as a function of calcination temperature; color chromaticity diagram (b); PL-PLE contour plots of a-YAB powders calcined at 470 1C (c) and 600 1C (d).

  (a) Material synthesisThe a-YAB powders with chemical composition close to the YAl 3 (BO 3 ) 4 compound were obtained by the sol-gel synthesis

Fig. 6

 6 Fig. 6 Time-gated phosphorescence contour plots of SG-500 (a) and SG-570 (b) a-YAB powders. The data in (a) represent the phosphorescence signals at 1 s after ceasing the excitation, whereas it is at 0.5 s in (b) as the afterglow intensity is poor for SG-570 due to a higher calcination temperature (ESI, † Fig. S4). The phosphorescence spectral profiles of a-YAB powders under 255 nm and 350 nm excitations are presented in (c).

Fig. 7

 7 Fig. 7 Allowed electronic transitions and energy levels of the first few prominent emitting centers in a-YAB powders.

Table 1

 1 Color rendering index, R a , correlated color temperature, CCT, luminous efficacy of radiation, LER and luminescence quantum efficiency, QY, of a-YAB powders under 365 nm LED excitation

	Samples	CRI (R a )	CCT (K)	LER (lm W À1 )	iQY (AE5%)
	SG-450	-	4100 000	181	32
	SG-470	82	22 600	215	31
	SG-500	85	9900	252	30
	SG-530	87	6450	261	34
	SG-570	86	5090	275	34
	SG-600	89	4370	264

Acknowledgements

This work was carried out under the ANR LuminoPhor-LED project (ANR-14-CE05-0033). The authors would like to thank Prof. Philippe Smet, UGent, for the discussion on phosphor properties.

Notes and references

Conflicts of interest

There are no conflicts to declare.


Vanni Petrolli

Magali Le Goff

Monika Tadrous

Kirsten Martens

Cédric Allier

Ondrej Mandula

Lionel Hervé

Silke Henkes

Rastko Sknepnek

Thomas Boudou

Giovanni Cappello
email: giovanni.cappello@univ-grenoble-alpes.fr

Martial Balland
email: martial.balland@univ-grenoble-alpes.fr

Confinement-induced transition between wave-like collective cell migration modes

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Supracellular organization plays a key role in establishing and maintaining structure, function and homeostasis in tissues. In the early stages of embryonic development, where features need to arise spontaneously from a homogeneous state, this organization closely follows morphogenic chemical patterns. Turing pioneered the idea that symmetry breaking results from linear instabilities in the reaction-diffusion dynamics of morphogens [1]. In the most general case, however, chemical reactions, osmotic pressures and mechanical forces all cooperate to determine the tissue-level organization. This is confirmed by an increasing number of recent studies, which indicate that cell proliferation, differentiation and motility are strongly impacted by the physical properties of the microenvironment [2][3][4][5][6]. Therefore, a full physical understanding of tissue mechanics and morphogenesis requires treating chemical and mechanical effects simultaneously. Several recent works reported that wave-like patterns of the local cell velocity spontaneously appear in colonies of epithelial cells [7,8]. Those velocity waves have also been observed in spreading epithelial sheets [9][10][11], regardless of cell proliferation [12], and are correlated to oscillations of the forces exerted by the cells on the substrate [13]. Such long wavelength patterns also appear in confined geometries where cell migration is limited to local cell rearrangements [14][15][START_REF] Duclos | Cell Migration[END_REF][START_REF] Deforet | [END_REF][18]. These waves are characterized by a wavelength λ and a period T , and show a surprisingly large spatial and temporal coherence. They can be modelled either at the particle level [START_REF] Deforet | [END_REF] or using continuum approaches [12,18], based on a coupling between cell motility and intercellular forces. In this Letter, we explore whether period and wavelength of collective wave excitations in epithelial cell monolayers are intrinsically encoded in the activity of the cell, or if they are affected by external constraints such as a specific set of boundary conditions. To achieve this, we analyzed the collective motion of epithelial cells confined to a quasione-dimensional channel. The experiments were accompanied by a series of numerical simulations, based on a self-propelled Voronoi model (SPV) [19][20][21], adapted to take into account the confining geometry. Our results show that tuning the length of the confining channel drives a phase transition between a state of global oscillations and a multinodal wave state. This transition is a consequence of the interplay between local cell active dynamics and global confinement. The effect is robust and does not require detailed knowledge of molecular processes but relies on a simple polarity-velocity alignment mechanism studied in the physics of dense active matter systems.

To confine cells to a quasi one-dimensional pattern, we prepared adherent stripes on soft polyacrylamide gels (E 40 kPa), as described previously [22] (outlined in Fig. 1a). Stripes of different length (100 to 2000 µm), but of the same width (40 µm), were patterned on the same substrate. Epithelial Madin-Darby Canine Kidney (MDCK) cells were then seeded on the patterned substrates with initial concentration of 2.5 ± 0.5×10 4 cells/cm 2 . The samples were washed with fresh medium 1h after seeding, then placed in the incubator (37 • C and 5% CO 2 ) until the end of the experiments. Cells were imaged in-situ unsing in-line holographic (defocus) microscope (see SI and Fig SI -1) [23] for 48 hours after confluence, gathering one image every 10 minutes (e.g., Fig. 1a-Phase). Cell velocities were computed with a custom-made Particle Image Velocimetry algorithm (PIV) with a final resolution of 20 min and 14 µm. To generate the kymograph, we cropped the videos in time to consider only confluent tissues, in an interval where the average absolute velocity was higher than 4 µm/h [24]. We then averaged the horizontal component of the speed along the transverse direction v x (x; t) = v (x, y; t) y . We removed low frequency drifts using a Gaussian highpass filter cropping 50% of the signal at 700 µm and 10 h. The kymograph in Fig. 1b represents the spatio-temporal evolution of the velocity field over 22 hours and over the whole stripe. Fig. 1c shows the instantaneous velocity profile, where a periodic oscillation of the speed appears. To quantify the period and the wavelength of these oscillations, we computed the autocorrelation function of the kymograph g(δx, δt) = v x (x, t) v x (x + δx, t + δt) x,t , displayed in Fig. 1d. We observe the appearance of an oscillating pattern in the autocorrelation function, both along the spatial and the temporal directions (respectively shown in panels e and f in Fig. 1). This pattern indicates the establishment of an extended multinodal standing wave, with the wavelength and period λ = 370 ± 30µm and T = 4.7 ± 0.7 h, respectively (errors represent the standard deviation, n = 59) (see histograms in Fig. 1e-f).

To obtain a detailed understanding of oscillations in tissues, we consider a computational framework based on a recently introduced self-propelled Voronoi model (SPV) [19][20][21]. The model used in this study is similar to that used in Ref. [21] to describe flocking transitions in confluent tissues, with an important modification that instead of using periodic boundary conditions, we imposed confinement through a repulsive rectangular wall of size (L X , L Y ), to reproduce the experiments geometry. Full details of the model and its implementation can be found in Ref. [20] (also see Supplemental Material for the parameters used). Briefly, the confluent cell monolayer is modeled as a two-dimensional network of Voronoi polygons covering the plane (Voronoi tessellation of all cell centre positions, see Fig. 2a). Each configuration of cells is described by the positions of cell centroids with its energy given by the commonly used Vertex Model [25], which depends on the area and perimeter of each cell. The parameters of the Vertex model include area and perimeter stiffness constants (K and Γ ) and target area and perimeter (A 0 and P 0 ). These parameters were chosen to describe a monolayer in a solid like regime (with a shape factor p 0 = P 0 / √ A 0 = 2.5) [19,26], to avoid shear flows induced by the boundaries. As in Refs. [19-we consider an overdamped dynamics, i.e., a force balance between frictional force with the substrate, selfpropulsion at a constant velocity v 0 along the direction of cell polarity, n i , and mechanical forces between the cells determined as a negative gradient with respect to cell position of the SPV model energy functional. The value of v 0 can be set to match the experimental observations (Fig. 2b-d), but does not affect the general oscillatory behavior. The dynamics of the cell polarity n i , described by the angle θ i with the x-axis of the laboratory reference frame (i.e., n i = (cos(θ i ), sin(θ i ))) is

∂θ i ∂t = 1 τ al sin(θ i -φ i ) + ν r i (t), (1)
with φ i being the angle between the velocity of cell i and the x-axis, and ν r i (t) being an orientational Gaussian noise. The angular dynamics is thus controlled by the interplay of rotational diffusion (kept constant in this study) and the polarity-velocity alignment with rate τ -1 al , with τ al being the time required by the cell to reorient its polarization in the direction of its velocity. This feedback mechanism leads to oscillations in confinement, where τ al plays the role of an effective inertia, and the oscillations are along the lowest-energy elastic modes of the material [START_REF] Henkes | Physical review[END_REF]. This feedback mechanism is also at the origin of flocks of in non-confined tissues [21].

Simulations of confined tissue layers show steady state oscillations akin to those observed in experiments (Fig. 2b) In the following, we study the dependence of these oscillations on the confining length L X and show that a feedback mechanism for alignment (through τ al ) is key to observe such mechanical waves in SPV. First we consider the case of long confining channels, where multinodal oscillations were observed in the experiments (Fig. 1). The simulation results displayed in Fig. 2c-d are obtained for a system with the same confining length (about 3 cells in y-direction) and aspect ratio as in the experiments in Fig. 1. We observe a pattern in the xcomponent of the velocity, v x , and using the same analysis tools as in Fig. 1, we extract the wavelength λ SP V ≈ 22 model length units and the period T SP V ≈ 8 model time units. Note that by approximately matching the timescale of the model to the experiments (through the cells velocity), one would get from these simulation data λ ≈ 300 µm and for the period T ≈ 2 hours. This indicates that this model is able to reproduce the features observed in the experiments, although some fine tuning of parameters is required for a quantitative match. Note that although the instantaneous velocity profiles and autocorrelation plots appear to be similar to the experiments, the full spatio-temporal dynamics of the model does not corresponds to standing wave oscillations. If the system size L X is decreased (keeping the value of τ al constant), the number of nodes also decreases up to a point where the system size can only accommodate a single spatial period of oscillation, and the system reaches a regime of global oscillation, where the direction of motion of all cells is coordinated. This transition is illustrated in Fig. 2e. The feedback timescale also plays an important role as no oscillations are observed if τ al is too large (i.e., the noise dominates over the coupling), and the critical length at which one observes multi-nodal oscillations increases with τ al . In the small system regime, the period of oscillations increases linearly with the system size as previously reported [START_REF] Deforet | [END_REF]18], and also with τ al (until oscillations eventually vanish for large values of τ al ), consistent with the role of the feedback mechanism as an effective inertia [START_REF] Henkes | Physical review[END_REF]. Therefore, the SPV model describes the transition controlled by the stripe length L c (τ al ) between global oscillations where all cells coordinate their motion to a regime where patches of cells coordinate their motion direction locally.

To verify this prediction, we varied the length L X of the stripe between 100 µm and 2000 µm (Fig. 3a), in order to tune the system across the critical length L c . In approximately 95% of experiments, in agreement with model predictions, we observed two types of behaviors: 1) A global movement of all cells that alternates between rightward and leftward motion (as seen from the autocorrelation function of the kymograph in Fig. 3b) and 2) The establishment of a multinodal standing wave with the antinodal cells moving back and forth, while cells in the nodes are being alternately compressed and dilated.

The incidence of the two behaviors strongly depends on L, with a transition for L λ. In the experiments with L < 200 µm, the global oscillation statistically dominated. In this case, the period scales linearly with the tissue size (Fig. 3d, blue area), while the wavelength is imposed by the confinement. In large structures (L > 500 µm), we only found multinodal waves, with the period and wavelength independent of L (Fig. 3d, red area). Fig. 3c quantifies the transition, with on average 39 tissues per point, obtained from three independent experiments. Our experiments confirmed the existence of a self-sustained oscillatory mode in the epithelial layers. Using the typical period and wavelength, we can define an effective velocity u φ = L/T 78 ± 13 µm/h, which is independent of the pattern size. Even for small patterns (L <500 µm), this velocity is preserved as the period scales linearly with the pattern length. We also note that u φ is approximately ten-fold larger than the average speed of individual cells within the epithelial layer (between 4 and 12 µm/h, depending on cell density [24,28]). Eventually, the spatial coherence of supracellular waves exceeds the largest pattern observable with our microscope.

Simulations using the SPV model show the emergence of sustained collective oscillations in confined monolayers. We identified two crucial conditions to produce these oscillations: 1) The existence of a delayed feedback between the cell velocity and self-propulsion direction that introduces a new timescale in the dynamics and 2) A very limited number of cellular rearrangement, at the limit of the solid-like regime. These ingredients allow the system to be described by linear elasticity, and for oscillations along the lowest energy elastic modes to dominate the dynamics [START_REF] Henkes | Physical review[END_REF]. One could thus envision tuning the oscillations by controlling cell-cell interactions through RAB5 or cadherin-mediated junctions, without affecting cells' individual mobility [21,29]. Contrarily to experiments, where multi-nodal standing waves are observed, the SPV model describes propagating oscillations. Several reasons may explain the difference. First, a standing wave is established only when the wavelength exactly matches the boundary conditions. Thus, a fine-tuning of the pattern length is needed in models, while the intrinsic variability between cells could make the real epithelium more adaptable to small variations of the confinement size. Second, a different choice of the coupling mechanism could also introduce a new timescale in the model and better describe standing waves in confined tissue. Two-dimensional SPV models are usually adapted to describe spatially extended monolayers, while the stringent confinement used experimentally makes the system be quasi-onedimensional. Moreover, whereas in the experiments the cells appear to be more elongated near the boundary, the constraint of maintaining a Delaunay triangulation (dual of the Voronoi tesselation) in the SPV prevents elongated cell configurations, and hence the tis- sue configuration in simulation do not realistically capture the shape of cells near the containment boundary. In conclusion, we demonstrate that the typical period and wavelength of epithelial tissue oscillations are intrinsically encoded in the cells, and are not adapted to external confinements. For this system, our SPV model predicts a transition between global oscillation and multi-nodal waves, the existence of which is confirmed experimentally for a pattern length L c 400 µm. From a biological perspective this transition could be significant. If in small systems all the cells behave similarly -the entire layer alternately moves back and forth -in large systems cells located either in the nodes or in the anti-nodes experience different mechanical stimuli and may undergo different fates, which can ultimately lead to supracellular patterning. The existence of an intrinsic wavelength λ also provides an intrinsic metric, likely encoded in the cell. It is interesting to note that λ roughly corresponds to the typical sizes of an embryo (e.g. in Drosophyla embryo, both length and circumference approaches 400-500 µm). Based on this, two important biological questions arise. Is this intrinsic metric used by the organism to measure distance inside a developing embryo? Does a collective long range excitation allow cells to probe their distant environment, in a timescale much shorter that that allowed by their own motility?

The authors would like to acknowledge I. Wang for the development of Particle Image Velocimetry and and P. Moreau for his technical support. We also thank K. Hennig, T. Andersen and C. Guilluy for valuable and supportive discussions. GC has been supported by the Institut National de la Santé et de la Recherche Médicale (Grant "Physique et Cancer" PC201407). T.B. and M.B. acknowledge financial support from the CNRS "Mission pour l'Interdisciplinarité" and the Center of Excellence of Multifunctional Architectured Materials "CEMAM" (n AN-10-LABX-44-01). M.B. acknowledges financial support from the ANR MechanoSwitch project, grant ANR-17-CE30-0032-01 of the French Agence Nationale de la Recherche. RS and SH the UK BBSRC, award numbers BB/N009789/1 and BB/N009150/1-2.

  FIG. 1. (a) Top: MDCK cells are seeded onto a polyacrylamide (PA) gel patterned with fibronectin stripes (width: 40 µm, length: 100-2000 µm).Center: Phase-contrast image of confluent tissue confined to a 1500 µm long stripe and the relative velocity field measured by PIV (bottom). Velocities pointing in the positive (negative) x-axis direction are shown in blue (red), in agreement with the arrows reported under the image. b) Kymograph representing the average horizontal velocity v (x; t) in time and and an example of velocity profile (c) along the dashed line. We removed low frequency drifts using a Gaussian high-pass filter. To quantify the periodicity of oscillations, we calculate the spatio-temporal autocorrelation of the kymograph (d) and measure peak spacing along the spatial (e) and temporal coordinates (f) (insets: distribution of peak periodicity for n=59 independent stripes). Images in panels (d) and (f) were smoothed for visualization purposes with a low-pass Gaussian filter (σx=15 µm, σt=30 min). Smoothing was not applied to panels (e), (g) and (h).

2 .

 2 Self-propelled Voronoi Model for collective oscillations in confluent tissues. (a) Example of tissue configuration obtained from the integration of the SPV model. Voronoi tessellation of the plane and centroid positions. (b) Velocity field of the centroids of the tesselation. Velocities pointing to the positive direction on the x-axis are represented in red and to the negative direction in blue. (c) Kymograph representing the average horizontal velocity (v (x; t)) over time and its profile (d) along the dotted line. (e) Phase diagram of oscillation patterns in the monolayer. Oscillations are observed for small enough values of the feedback timescale (τ al < τ c al ≈ 17 model time units). For τ al < τ al , multi-nodal oscillations are observed for large systems LY > Lc(τ al ) whereas small systems oscillate as a block (global oscillations). The dots delimiting the two phases indicate the value of Lc(τ al ) and are determined by the number of nodes in the spatial velocity autocorrelation functions (0 in the case of global oscillations) and the full line is a power law fit. Insets show the autocorrelation of the typical kymographs of x-velocity in the corresponding phases.

FIG. 3 .

 3 FIG.3. Dependence of oscillatory behaviour on the stripe length. a) The velocity field superimposed on phase contrast images for short stripes of length 200 µm and 300 µm displays global oscillations, generating a characteristic twodimensional autocorrelation (right). Longer lines (500 µm and 1000 µm) display multinodal oscillations (b), which give rise a different pattern in the autocorrelation image. Velocities pointing in the positive x-axis direction are represented in blue, those pointing in the negative x-axis direction are represented in red, in agreement with the arrows reported in the schemes under each image. For each length, we display the frequency of each phenotype (c) and the characteristic time and space periodicity (d) calculated. Bars represent the standard error of the mean.


Monika E Dolega

Benjamin Brunel

Magali Le Goff

Magdalena Greda

Claude Verdier

Jean-François Joanny

Pierre Recho
email: pierre.recho@univ-grenoble-alpes.fr

Giovanni Cappello
email: giovanni.cappello@univ-grenoble-alpes.fr

Affiliations

Extracellular Matrix acts as pressure detector in biological tissues Authors

 , with a compliance that has to match the perturbation amplitude.

 . We support this finding with a theoretical model of an actively pre-stressed composite material, made of incompressible and impermeable cells and a poroelastic interstitial space. This description correctly predicts the emergent bulk modulus of the aggregate as well as the hydrodynamic diffusion coefficient of the percolating interstitial fluid under compression. We further show that, on a longer timescale, the extracellular matrix serves as a sensor that regulates cell proliferation and migration in a 3D environment through its permanent deformation and dehydration following the global compression.

Introduction

There are several evidences that the mechanical environment of a cell plays a fundamental role on cell fate, both in terms of proliferation and differentiation. Discher et al. [START_REF] Engler | Matrix Elasticity Directs Stem Cell Lineage Specification[END_REF] observed that stem cells differentiate into different cell types, when cultured on polyacrylamide gels of different stiffness covered with collagen. Trappmann et al. [START_REF] Trappmann | Extracellular-matrix tethering regulates stem-cell fate[END_REF] ascribe such a differentiation to the collagen/polyacrylamide cross-link density, which accidentally correlates with the gel stiffness. Nevertheless, both studies point to the crucial role played by the mechanical interaction between the cell and its microenvironment. This becomes particularly relevant in an tissue structure, where the cells are in contact with an ExtraCellular Matrix (ECM) which is essential to the structure and function of the tissue as cells respond to chemical and mechanical signals provided by the surrounding ECM, but also constantly remodel and/or produce it [START_REF] Bonnans | Remodelling the extracellular matrix in development and disease[END_REF] . ECM biochemical composition and 3D organization varies from thin layers of basement membrane, in epithelial tissues [START_REF] Leblond | Structure, composition, and assembly of basement membrane[END_REF] , to more abundant layers, in between cells in tumors [START_REF] Levental | Matrix Crosslinking Forces Tumor Progression by Enhancing Integrin Signaling[END_REF][START_REF] Jodele | Modifying the soil to affect the seed: role of stromal-derived matrix metalloproteinases in cancer progression[END_REF] . Such abundant and crosslinked ECM leads to tissue stiffening, which in turn increases the risk for malignancy by enhanced integrin signaling [START_REF] Levental | Matrix Crosslinking Forces Tumor Progression by Enhancing Integrin Signaling[END_REF] . Yet, how the presence of the ECM in multicellular structures determines their mechanical behavior remains to be determined and could be a key factor in some cases. This is exemplified by the fact that the compressibility modulus of individual cells [START_REF] Monnier | Effect of an osmotic stress on multicellular aggregates[END_REF] is several orders of magnitude higher than the one of multicellular aggregates of similar cells [START_REF] Netti | Role of extracellular matrix assembly in interstitial transport in solid tumors[END_REF][START_REF] Leroux | Elastography of multicellular aggregates submitted to osmomechanical stress[END_REF]

Results

To quantify the ECM interstitial volume fraction 𝑛 𝑚 in MCS, we performed negative imaging by supplementing the culture medium with sulphorhodamine B, a hydrophilic fluorophore that permeates the luminal space (Figure 1b). This approach allows for live confocal imaging and avoids the potential artefacts due to fixation and cryosectioning. Taking into account the confocal resolution (270 nm in our experimental conditions),

we measured an average thickness of the intercellular space of 0.9±0.1 µm. Considering that CT26 cells have a typical diameter of about 20µm, we estimate an interstitial volume fraction of 𝑛 𝑚 = 14±5% (Figure 1c and SI-1).

To understand the role of ECM in perceiving and transmitting mechanical perturbations we established a novel method that allows compressing either globally our model MCS, or selectively only the cells, leaving the ECM unstrained (see sketch in Figure 2a).

Dextran molecules with small hydrodynamic radii (r < 6 nm, MW < 70 kDa) penetrate the ECM, without entering the cells [START_REF] Cadart | Size control in mammalian cells involves modulation of both growth rate and cell cycle duration[END_REF] and thus build up a hydrostatic pressure across their membranes that selectively compresses them. Instead, a global compression is obtained by exposing MCS to an osmotic stress exerted by big Dextran molecules (r > 15 nm, MW > 200 kDa) that are excluded from both cells and ECM (Figure 2b-left) [START_REF] Montel | Stress Clamp Experiments on Multicellular Tumor Spheroids[END_REF][START_REF] Monnier | Effect of an osmotic stress on multicellular aggregates[END_REF] .

First, we assessed the long-term consequences of such global versus selective osmotic compression using three parameters: cellular proliferation, cellular motility and MCS roughness. We observed that whereas cells in control MCS (0 kPa) present a rather uniform proliferation pattern, a global compression of MCS (big Dextran) stops cell division in the core and alters the overall MCS growth (Figure 2f), as previously reported [START_REF] Helmlinger | Solid stress inhibits the growth of multicellular tumor spheroids[END_REF][START_REF] Alessandri | Cellular capsules as a tool for multicellular spheroid production and for investigating the mechanics of tumor progression in vitro[END_REF][START_REF] Montel | Isotropic stress reduces cell proliferation in tumor spheroids[END_REF] . Such altered cellular proliferation under 5 kPa resulted in a three times slower MCS growth as compared to control (Figure 2g). In parallel, cell motility within the MCS at 5 kPa (big Dextran) decreased by 50% (Figure 2i and SI-2).

We also observed that global compression induces MCS smoothening (Figure 2f bottom and Figure 2h), which we assume is related to an increased cell cohesion. Strikingly, all three parameters were not altered when MCS were exposed to an equivalent pressure (small Dextran) applied selectively to cells residing in a native ECM, suggesting the major role of ECM in transducing changes in mechanical state of the environment to cells. However, the ECM properties under compression are difficult to characterize inside the MCS directly.

To probe the rheological properties of the ECM upon osmotic compression, we used Matrigel (MG), an ECM proxy secreted by EHS mouse sarcoma [START_REF] Kleinman | Isolation and characterization of type IV procollagen, laminin, and heparan sulfate proteoglycan from the EHS sarcoma[END_REF] .

Consistently with native ECM, large Dextran

Discussion

Several hypotheses can be put forward as to how cells transduce such residual stress into a biological response. First, MG compression is associated with its stiffening, as the Young's modulus increases 6-folds (Figure 2e). This would be coherent with the evidence that the substrate stiffness influences the cellular fate [START_REF] Engler | Matrix Elasticity Directs Stem Cell Lineage Specification[END_REF] . Compression also results in a reduction of MG porosity by ~25% (Figure 2c and SI-1), with the consequent reduction of diffusion of oxygen, nutrients and metabolites as well as insoluble cues present in the ECM (e.g. GPB, cRGD) [START_REF] Wrighton | Signals from the surface modulate differentiation of human pluripotent stem cells through glycosaminoglycans and integrins[END_REF] . However, for single cells seeded in MG subjected to a global compression, the invasive morphotype was not rescued with MG of higher densities corresponding to the compressed state (Figure 3c), suggesting that the presence of the stress in the external environment directly induces the morphotype change rather than its stiffness or porosity.

Second, ECM compression may restore cell-cell contact inhibition of proliferation 8 and locomotion [START_REF] Garcia | Physics of active jamming during collective cellular motion in a monolayer[END_REF] through maturation of cell-cell adhesion. Indeed, as

shown by phase contrast images in Figure 2f and quantified in Figure 2i, the whole MCS compression smoothened its surface, which we interpret as an increase of cell-cell adhesion. Conversely smoothing did not occur with small osmolites. The fact that mechanical compression arrests cell motility is also in agreement with a simple model of the cell cytoskeleton self-organization [START_REF] Putelat | Mechanical stress as a regulator of cell motility[END_REF] .

Interestingly, a recent study points to the role of the ECM swelling in tumor growth [START_REF] Voutouri | Hyaluronan-Derived Swelling of Solid Tumors, the Contribution of Collagen and Cancer Cells, and Implications for Cancer Therapy[END_REF] . Because the ECM of tumors contains negatively charged macromolecules, counter-ions permeate the MCS to establish electroneutrality and create an osmotic pressure difference between the MCS inside and outside. This pressure difference is equilibrated by an interstitial hydrostatic pore pressure within the MCS that swells the ECM. It is then conceivable that such deformation, opposite to the ECM compression that we measure during osmotic compression, sustains the MCS growth. Such an effect can be included in our theoretical model by accounting for the ionic species in solution in the MCS pores [START_REF] Xue | A nonlinear poroelastic theory of solid tumors with glycosaminoglycan swelling[END_REF] , however we neglected kDa; Sigma-Aldrich, St. Louis, MO) was added to the culture medium to exert mechanical stress, as previously described [START_REF] Monnier | Effect of an osmotic stress on multicellular aggregates[END_REF] , at a concentration of 55 g/L to exert 5 kPa.

Materials and Methods

Cell

To follow spheroid growth over the time, phase contrast images were taken daily starting from "day 0" before addition of dextrans and after 30 minutes.

Matrigel beads preparation

Matrigel beads were prepared using vortex method.

Oil phase of HFE-7500/PFPE-PEG (1.5 %w/v) was cooled down to 4°C. For 400 µL of oil, 100µL of Matrigel was added. Solution was vortexed at full speed for 20 seconds and subsequently kept at 37'C for 20 minutes for polymerization. Beads were transferred to PBS phase by washing out the surfactant phase with pure HFE-7500 oil.

Matrigel beads compression

To compress polimerized Matrigel beads, PBS was enriched with either 2MDa or 1kDa dextran at the concentration to exert a pressure of 5 kPa. Images were taken just before, and then 45 minutes after, dextran was added. Volume decrease was measured for 10 different beads.

Matrigel

Infusion of dextrans into Matrigel beads

18 hours before observation, Fluorescent dextrans were added to PBS containing Matrigel beads (in control or compressed as previously described) to the final concentration of 500 nM for 500kDa dextran, and 5 µM for 40kDa and 70kDa dextran.

Statistical analysis

Student

Supplementary Materials

See attached file.

Figures and Tables

  . In this article, we use composite MultiCellular Spheroids (MCS) composed of murine colon cancer cells (CT26 WT) to design an experimental technique to either selectively compress the single cells in the spheroid or globally the whole MCS. Next, we propose a mechanical model that quantitatively predicts the emergent compressibility and permeability of the MCS in which impermeable and incompressible cells coexist with a soft poroelastic ECM (imaged by immunolabeling of fibronectin in Figure 1a) permeated by interstitial fluid. In parallel, we show that this peculiar rheology has a strong impact on cell proliferation and migration in MCS under compression.

  molecules were also excluded from microbeads made of MG (Figure 2b-right and SI-1) suggesting an equivalent effective permeability. Upon 5kPa compression, MG microbeads shrunk by approximately 64±5% (Figure 2c-d), which indicates that even small mechanical perturbations in vivo are capable to impose important changes to the residing ECM. As measured by atomic force microscopy (AFM), the Young's modulus of compressed MG gels is significantly larger (610±70 Pa) as compared to control (103±4 Pa), suggesting a strain stiffening behavior of MG (Figure 2e and SI-3). We observed a significantly smaller compression (10±4%) with small Dextran, which enters the MG 27 . Using this experimental system, we applied compression to individual cells seeded in MG to reproduce a 3D-like environment with no neighbors and confirm our results obtained for whole MCS in order to rule out collective effects stemming from potential cell-cell contacts. Few days after insemination in MG, the control cells (0 kPa) presented long protrusions, reminiscent of the typical morphology of invasive malignant cells. A global compression of MG with large osmolites resulted in strikingly different nonprotrusive morphotype, while a selective compression of the same magnitude applied via small osmolites was insufficient to alter the morphotype (Figure 3a and 3b). Our results therefore suggest that the ECM plays the role of a mechanical proxy through which the global MCS compression impacts the cells fate. To understand the mechanism behind such mecanosensitivity, in line with experimental observations shown on Figure 1 and Figure 2, we modelled the short time scale response (i.e. associated to water percolation) of the MCS to a global osmotic compression using a composite description, where cells are impermeable and incompressible and the ECM is poroelastic (Figure 4b). A pre-stress, due to the cell division and loss within the MCS, was also taken into account consistently with 28 and assumed to be fixed at the timescale of the compression. We obtained (see details in SI-4) the following formulas for the effective volumetric bulk modulus and the hydrodynamic diffusion coefficient associated to water percolation during compression: the above relations, 𝜇 is the viscosity of the extracellular fluid, 𝐾 𝑚 is the bulk modulus of the ECM, 𝜅 𝑚 the permeability of the ECM, 𝑝 𝑎 is the active stress in the MCS and 𝐺 its shear modulus. Rigorously, such formulas are derived in a linear theory assuming that 𝑝 𝑎 𝐺 ⁄ is small. However, we shall still use them for moderate values of 𝑝 𝑎 𝐺 ⁄ ≈ -0.7 as estimated by fitting the theoretically predicted profiles to the experimentally measured local stress and strain 28,29 within the MCS (details in see SI-4.3). As we measure by AFM a shear modulus of 𝐺 = 280 ± 10 Pa (Figure d and SI-3.2), we thus estimate an active pressure of the order of 𝑝 𝑎 =-200 Pa. Both 𝐾 𝑚 = 1.2 ± 0.2 kPa (Figure c and SI-4.4) and 𝐷 𝑚 = 30 ± 3 𝜇m 2 /s (Figure e and SI-4.4) were directly obtained from the measurement of the compressibility and relaxation time of MG beads under osmotic pressure. Together with the measured value of 𝑛 𝑚 , formulas (1) have no adjustable parameter and predict a value of 𝐾 of few tens of kPa and 𝐷 of few tens of 𝜇m 2 /s. Such values are consistent with direct measurements of 𝐾= 29±3kPa (Figure c and SI-4.3) and 𝐷 = 44±4 𝜇m 2 /s (Figure e and SI-4.3) based on global compression experiments confirming the model applicability. In this framework, the ECM plays the role of a sensor transmitting the small compressive stress to the cells in the aggregate. Indeed, a global compression increases the hydrostatic pressure in the inter-cellular space, draining the water out and straining the ECM (Figure 4f), which in turn mechanically compresses the cells in a permanent way as such stress would be relaxed only with complete remodeling of the ECM. In contrast, cells rapidly (1-10min) offset a small osmotic stimulus -5 kPa (this corresponds to 1 mM NaCl, which is a very small variation compared the physiological values of concentration of [Na],[Cl]≈150 mM in the extra-cellular space) through a regulatory response, which activates ion pumps 30 to equilibrate the internal/external osmolarity unbalance. As a consequence, an omosotic shock selectively compressing the cells has a fundamentally different biological signature on MCS than a global osmotic shock, which exposes the cells to a permanent mechanical pressure. The model also indicates that water percolation through the MCS pores effectively follows a diffusive process. Thus, the timescale at which the incompressible/compressible transition occurs scales quadratically with the typical size. In vision of that, tissues are incompressible because of the presence of water but, with the measured value of 𝐷, an individual cell (typically 20 μm) 'feels' an incompressible environment only for few minutes if it swells or moves. As a result, the slow processes of cell division and motion take place in an effectively compressible environment, but where some residual stress stemming from the ECM deformation remains until the ECM is completely remodeled.

  it given that high (200 mM NaCl) modifications of counter ions concentrations are necessary to record a non-negligible pore pressure increase while we imposed much smaller osmotic perturbations (5 kPa corresponds to 1 mM NaCl). The viability of cells in such osmotic environment would need to be carefully scrutinized.To conclude, our experiments are compatible with the model proposed here, where ECM mechanical properties along with volume exclusion due to incompressible cells determines the MCS final volume loss and its dynamics upon a gentle global compression. In this framework, the ECM is a sensor through which the global MCS compression is transformed into a permanent biochemical signal impacting proliferation (Figure4a). While we have shown that the MCS rheology and, in particular, the state of compression of the ECM affect cellular fate within the MCS, cellular turnover will in turn modify the active pre-stress leading to an emergent hydrodynamic diffusion and mechanical stress within the MCS in the long timescale. We therefore anticipate that mechanical theories aiming at capturing such state[START_REF] Ambrosi | Solid Tumors Are Poroelastic Solids with a Chemomechanical Feedback on Growth[END_REF] should further account for the presence of the ECM.

  culture, MCSs formation, and growth under mechanical stress CT26 (mouse colon adenocarcinoma cells, ATCC CRL-2638; American Type Culture Collection were cultured under 37°C, 5% CO2 in DMEM supplemented with 10% calf serum and 1% antibiotic/antimycotic (culture medium). Spheroid were prepared on agarose cushion in 96 well plates at the concentration of 500 cell/well and centrifuged initially for 5 minutes at 800rpm to accelerate aggregation. After 2 days, Dextran (molecular mass 1, 10, 40, 70, 100, 200, 500 and 2000

  Spheroid were kept under constant pressure over observation period. Images are analysed using Imagej plugin (doi.org/10.1371/journal.pone.0103817). Each experiment was repeated 3 times, with 32 individual spheroids per condition. Spheroids were fixed with 5% formalin (Sigma Aldrich, HT501128) in PBS for 30 min and washed once with PBS. For cryopreservation spheroids were exposed to sucrose at 10% (w/v) for 1 hour, 20% (w/v) for 1 hour and 30% (w/v) overnight at 4°C. Subsequently spheroids were transferred to a plastic reservoir and covered with Tisse TEK OCT (Sakura) in an isopropanol/dry ice bath. Solidified samples were brought to the cryotome (Leica CM3000) and sectioned into 15µm slices. Cut layers were deposited onto poly-L-lysine coated glass slides (Sigma) and the region of interest was delineated with DAKO pen. Samples were stored at -20°C prior immunolabelling. For fibronectin and Ki67 staining samples were permeabilized with Triton X 0.5% in TBS (Sigma T8787) for 15 minutes at RT. Nonspecific sites were blocked with 3% BSA (Bovine serum Albumin) for 1 hour. Then, samples were incubated with first antibody (Fibronectin, Sigma F7387, 1/200 and Ki67; Millipore ab9260, 1/500) overnight at 4°C. Subsequently samples were thoroughly washed with TBS three times, for 15 minutes each. A second fluorescent antibody (goat anti-mouse Cy3, Invitrogen; 1/1000) was incubated for 40 minutes along with phalloidin (1/500, Alexa Fluor 488, Thermo Fisher Scientific). After extensive washing with TBS (four washes of 15 minutes) glass cover slides were mounted on the glass slides with a Progold mounting medium overnight (Life Technologies P36965) and stored at 4°C before imaging.

  porosity 18 hours before observation, fluorescent dextrans were added to PBS containing Matrigel beads (in control or compressed as previously described) to the final concentration of 500 nM for 500kDa dextran, and 5 µM for 40kDa and 70kDa dextran. Images were taken 30µm above the glass surface by the border of the Matrigel bead. Infusion efficiency was quantified by measuring the normalized Intensity (ImageJ) in at least 9 small regions within Matrigel (I gel ) and outside in the solution (I sol ). At least 8 images were taken per condition.

Figure 1 :

 1 Figure 1: Extracellular Matrix volume fraction in multicellular spheroids. (a) Immuno-fluorescent staining of fibronectin and (b) intercellular space observed by confocal microscopy. The intercellular space is filled by sulphorhodamin-B, a hydrophilic colorant that does not penetrate the plasma membrane. (c) Taking into account the confocal resolution, we measure that the intercellular space has an average thickness 930±50nm. We estimate that the extracellular space represents 14%±5% of the total volume of a whole multicellular spheroid.

Figure 2 :

 2 Figure 2 : Extracellular Matrix selective compression in MCS. (a) Selective compression method. Small Dextran molecules (blue) with MW < 70 kDa and a hydrodynamic radius smaller than 6 nm penetrate the ECM and exert a selective compression on the cells, but have no effect on the ECM. Conversely, big Dextran molecules (green) with MW > 200 kDa (R>15 nm) exert a global compressive stress on the whole MCS. (b) Left: Confocal sections of MCS dipped in a solution containing either small or big fluorescent Dextran molecules. The section shows that ECM is selectively permeable to small Dextran, with a cutoff radius of the order of 10 nm. Right: MG beads are also permeable to small Dextran, with a cutoff radius similar to that of ECM. (c) Compression of MG beads under 5 kPa, triggered by small and large Dextran molecules. (d) Beads lose 60% of their initial volume when compressed using big Dextran, and only 10% with small Dextran. Error bars represent the standard deviation of the mean. N=10. (e) AFM measurement of MG stiffness indicates that, under compression, MG stiffens by 6-fold. N=9 per each condition. (f) Above: Proliferating cells revealed by immunostaining of KI67 with no pressure, under overall compression of 5 kPa and under selective compression of the cells and the consequence on MCS growth speed (panel g), N > 36 MCS per each condition. Three independent experiments. Below: Effect of compression on the MCS roughness, defined as 𝑃𝑒𝑟𝑖𝑚𝑒𝑡𝑒𝑟/ 2𝜋√𝐴𝑟𝑒𝑎 and quantification in panel (h), N > 36 MCS per each condition. Three independent experiments. (i) Cell migration speed (see SI) within MCS also significantly depends on ECM compression. N= 5 independent experiments per condition. All results are displayed as bar graphs ± SD. Biological experiments were repeated at least three times.

Figure 3 :

 3 Figure 3: Compression of individual cells in MG. (a) MG compression promotes the growth of sphere-like clusters, while cells form stretched structures without pressure, or with a selective pressure but no MG compression. Above: Maximal projection of 50 µm Z-stack, actin staining; Below: analysis of shape anisotropy as a measure of circularity (C); high protrusive phenotype in magenta (C≤0.8) and cluster-like phenotype in cyan (C> 0.8) (b) Cumulative distribution of clusters anisotropy in control and under global and selective compression. (c) MG initial concentration has little effect on single cell/cluster anisotropy.

Figure 4 :

 4 Figure 4: Mechanical response of MCS and MG to mechanical stresses. (a) Schematic description of the experiment. The osmotic compression (magenta) is applied after 2 days of unconstrained growth (black curve). The MCS volume decrease typically lasts after few minutes, between marks (1) and (2). (b) The MCS is modeled as a composite material made of incompressible cells (bulk modulus K c ~ 1 MPa) surrounded by a porous ECM. ECM is characterized by its bulk and shear moduli (K m and G m ), its volume fraction n m and the porosity κ m . Marks (1) and (2) correspond respectively to the conformations before and after compression (see panel a) (c) The bulk modulus of MG and MCS (K m = 1.2±0.2 kPa and K = 29±3 kPa) and (d) the corresponding shear moduli (G m = 100±50 Pa and G = 280±10 Pa) measured by AFM as described in SI. (e) Diffusion coefficient of water through a bead of pure MG (D m = (3.0±0.2)•10 -11 m 2 /s) and through the composite MCS (D = (4.4±0.4)•10 -11 m 2 /s ), under compression. (f) Propagation of three different stresses in a composite poroelastic MCS. Left: Big Dextran exert the stress on cells through ECM compression (ε m > 0), which cannot be actively screened by cells over a short timescale. Proliferation decreases significantly. Middle: Small Dextran molecules exert the stress directly on the cells, but with no ECM deformation (ε m = 0), and has little impact on proliferation. Indeed the stress is suppressed over a long timescale by an active pumping of ions inside the cell to re-equilibrate the osmotic pressure. Right: ECM hydration observed in in-vivo tumors induces ECM distension (ε m < 0) and an overall MCS volume increase which promotes tumor growth.

  's t-test (unpaired, two tailed, equal variances)

	was used to calculate statistical significance as
	appropriate by using PRISM version 7 (graphpad
	Software). Statistical significance was given by *,
	P<0.05; **, P<0.01; ***, P<0.001; ****, P<0.0001.

Acknowledgments

We warmly thank J. Prost and F. Jülicher for drawing our attention to the potential impact of the poroelasticity on the rheology of multicellular aggregates.

Funding: . This work was supported by the Agence Nationale pour la Recherche (Grant ANR-13-BSV5-0008-01), by the Institut National de la Santé et de la Recherche Médicale (Grant « Physique et Cancer » PC201407) and by the Centre National de la Recherche Scientifique (Grant MechanoBio 2018). This work has been partially supported by the LabeX Tec 21 (Investissements d'Avenir: grant agreement No. ANR-11-LABX-0030).