Phototoxicity induced in living HeLa cells by focused femtosecond laser pulses: a data-driven approach

: Nonlinear optical microscopy is a powerful label-free imaging technology, providing biochemical and structural information in living cells and tissues. A possible drawback is photodamage induced by high-power ultrashort laser pulses. Here we present an experimental study on thousands of HeLa cells, to characterize the damage induced by focused femtosecond near-infrared laser pulses as a function of laser power, scanning speed and exposure time, in both wide-field and point-scanning illumination configurations. Our data-driven approach offers an interpretation of the underlying damage mechanisms and provides a predictive model that estimates its probability and extension and a safety limit for the working conditions in nonlinear optical microscopy. In particular, we demonstrate that cells can withstand high temperatures for a short amount of time, while they die if exposed for longer times to mild temperatures. It is thus better to illuminate the samples with high irradiances: thanks to the nonlinear imaging mechanism, much stronger signals will be generated, enabling fast imaging and thus avoiding sample photodamage. Here we present an experimental study of photodamage induced by focused ultrashort (130-fs) laser pulses on HeLa cells, under two different illumination configurations and varying laser power, scanning speed and exposure time. The results are presented following a data-driven approach, offering an interpretation of the underlying damage mechanism and providing a predictive model for biologically safe imaging in the range of parameters explored. Analyzing the dependence of the damage on the applied laser power, we mainly obtained photothermal damage, and we confirmed this result also with thermodynamic considerations. From our experiments, we calculated the probability of damaging the sample as a function of the variables examined, and we identified a safety that microscopy 5% damage


Introduction
Since its introduction [1], nonlinear optical microscopy has established itself as the tool of choice in many scientific fields, working in different modalities such as two-photon excited fluorescence, second/third harmonic generation and coherent Raman scattering. Also thanks to new laser sources that have become more and more compact, turn-key, stable, powerful and that generate ultrashort pulses with multi-wavelength outputs at high repetition rates, nonlinear optical microscopes are nowadays capable of directly imaging tissues and probing dynamic cellular processes in a label-free manner, reaching a remarkable spatial and temporal resolution [2,3]. The downside of these technological advances is the high-intensity illumination required, that can lead to phototoxicity. The general term phototoxicity includes all physical and chemical damaging processes that may arise during the interaction between light and biological organisms [4][5][6] and represents a serious limitation for in vivo biomedical imaging [7][8][9].
Phototoxic mechanisms can be induced by linear or nonlinear optical processes [10]. In the following, we present a brief overview of photodamage processes, starting from the linear mechanisms triggered at low light intensities and then describing the more destructive nonlinear mechanisms arising at higher and higher intensities. In laser physics, it is common to use the word intensity of radiant energy as the power transferred per unit area (W/m 2 ). In this manuscript, adhering to the language used in photometry, we will use the term irradiance to refer to the average laser intensity. Linear photodamage is mainly caused by sample heating due to laser irradiation [11], namely photothermal damage, but it may also be associated with photochemical damage [12]. This occurs in particular in pigmented molecules and pigment-rich tissues, such as skin [13][14][15], which may directly absorb light via single-photon transitions, especially when using visible or UV instead of near-infrared (NIR) radiation, even in the continuous-wave regime. Nonlinear photodamage mechanisms, induced by e.g. multiphoton absorption [16,17], arise using higher light intensities. The most common one is the photochemical damage correlated with the formation of reactive oxygen species (ROS) [18][19][20] and free radicals [21,22]. The formation of these toxic species is induced by nonlinear excitation of cellular and tissue absorbers, such as porphyrins, flavins and coenzymes as NAD(P)H. Under physiological conditions, cells have different mechanisms to cope with ROS production. However, when excess laser irradiation triggers ROS overproduction, not only cell components, such as proteins, lipids and DNA, begin to oxidize, but also redox homeostasis and cell cycle are disrupted [23]. Besides ROS production, cells may undergo direct DNA photodamage [22,24,25] and, if the irradiance is high enough, namely above ∼ 6 · 10 12 W/cm 2 for a NIR femtosecond beam, molecular ionization can lead to low-density plasma formation [21,26]. Lastly, when the plasma is no longer confined, optical breakdown occurs with immediate mechanical destructive effects for cells and tissues [27][28][29][30].
If the definition of phototoxic mechanisms is straightforward, quantifying phototoxicity during microscopy is a completely different ballgame [31]. In fact, photoxicity depends on many different factors, including sample type, specific sample sensitivity induced by experimental protocols and illumination parameters, both related to the experimental setup (continuous or pulsed source, illumination wavelength, irradiance, etc.) and to illumination protocols (imaging time, stationary or scanned beam, number of 3D planes measured, etc.) [29,[32][33][34][35]. Assessing the short-and long-term detrimental effects of laser radiation on biological specimens is a non-trivial problem and various methods have been elaborated. Among these, we mention measuring photobleaching and irradiation-induced fluorescence [9,36] and many in vitro mechanism-specific tests to target ROS and free radicals production, DNA strands, changes in methabolic activity and in cell morphology, stress-related proteins, cloning efficiency and level of apoptosis [37][38][39][40][41]. Despite the large number of studies conducted in this area [4,5,31,42,43], there are still many things to understand in order to formulate a general theory of photodamage, due to the enormous variety of parameters and the inability to explore them all.
For reasons of consistency with our work, in the remainder of this introduction we will focus on studies on label-free cultured cells. Among the most relevant results, Nan and coworkers [44] presented the possibility to perform nonperturbative Coherent Anti-Stokes Raman Scattering (CARS) microscopy. They used two pulsed laser beams, pump and Stokes, respectively at 711 nm and 892 nm, with a pulse duration of 2 ps, and they explored the occurrence of linear and nonlinear damage upon variation of pulse energy and average power. They characterized linear damage using low pulse energy (<1 nJ) and high average power, finding the linear damage threshold for average power between 9 mW and 12 mW. On the other hand, for nonlinear damage, they applied high pulse energy (>1 nJ) and average power below 5 mW. After 300 s of scan (100 frames), for 1 MHz of repetition rate they set nonlinear damage threshold for pump pulse energy between 2 nJ and 3 nJ. Repeating the same experiment for 2 MHz of repetition rate, they set nonlinear damage threshold for pump pulse energy between 1.5 nJ and 2.3 nJ. In both cases, the pulse energy of the Stokes was kept constant at 0.5 nJ. Eventually, they identified 1.5 mW (8 MHz, 0.2 nJ) for the pump beam and 0.5 mW (8 MHz, 0.06 nJ) for the Stokes beam as safety condition for nonperturbative CARS microscopy imaging. Also Fu and colleagues characterized photodamage in CARS microscopy [29], with the aim of producing a guideline to obtain high quality CARS images of biological samples avoiding photodamage. They analyzed the dependence of photodamage rate on excitation power, excitation wavelength, repetition rate. They also characterized the so called 'coherent Raman induced' photodamage, i.e. the photodamage caused by vibrational absorption in the sample in the presence of Raman resonance. Their results on living cells, employing a pump beam at 704 nm and a Stokes beam at 808 nm, show that for a total average power of 9.24 mW (peak power of 373 W and 101 W for pump and Stokes respectively) cells reported severe damage after 1 minute of exposure. Decreasing the total average power to 2.3 mW (peak power of 92.3 W and 24 W for pump and Stokes respectively), no damage was reported, even after 10 minutes of exposure. In their work on wide-field CARS microscopy, Silve et al. [45] proved the possibility to acquire CARS images of living cells in suspension. For their experiments, they used a pump beam at 625 nm and a Stokes beam at 770 nm, with pulse duration of 4 ns and 3 ns, respectively and 10 Hz of repetition rate. In their configuration, they demonstrated that no damage occurred after exposing the cells to 3000 pulses, with pump and Stokes pulse energy of 1.2 mJ and 0.53 mJ, respectively. In their work on stimulated Raman scattering (SRS) microscopy [46], Freudiger and coworkers accurately optimized pulse parameters for biological samples. Changing pulse duration, they experimentally identified three average power damage thresholds: 25 mW at 180 fs pulses, 80 mW at 1 ps pulses and 280 mW at 6 ps pulses. From these data, they extracted an empirical equation to evaluate the maximum permissible total average power to avoid photodamage, considering a pump beam at 817 nm and a Stokes beam at 1064, and a 1.2 NA objective.
Finally, Waldchen and colleagues [47] presented an extensive work on photodamage of living cells, analyzing three different cell lines and designing ad-hoc experiments to investigate the influence of irradiance, wavelength, light-dose, fluorescence labeling, temperature and reducing agent on the viability of labeled and unlabeled cells. Their experiments clearly showed how shorter irradiation wavelengths cause more damage: exposing unstained cultured cells to cw irradiation for 240 s at the same irradiance, 0.2 kW/cm 2 , the survival rate went from 0% to 100% changing the laser wavelength from 488 nm to 514 nm. At 640 nm, the irradiance can go up to 6 kW/cm 2 , without causing any damage to the irradiated sample. In addition, they proved the insufficiency of the light dose alone to estimate the degree of photodamage.
Here we present an experimental study of photodamage induced by focused ultrashort (130-fs) laser pulses on HeLa cells, under two different illumination configurations and varying laser power, scanning speed and exposure time. The results are presented following a data-driven approach, offering an interpretation of the underlying damage mechanism and providing a predictive model for biologically safe imaging in the range of parameters explored. Analyzing the dependence of the damage on the applied laser power, we mainly obtained photothermal damage, and we confirmed this result also with thermodynamic considerations. From our experiments, we calculated the probability of damaging the sample as a function of the variables examined, and we identified a safety window that can be used for microscopy applications below the 5% damage probability.

Cell culture, seeding and staining
For our experiments we used HeLa cells, a human cell line derived from cervical cancer cells, which are extensively employed in in-vitro phototoxicity experiments, and in particular in the field of laser microscopy, thanks to their availability, ease of cell manipulation and reproducibility of the results [47,48]. HeLa cells are cultured in Dulbecco's modified Eagle's medium (DMEM) supplemented with 10% fetal bovine serum, 1% L-glutamine (2 mM), penicillin (10 units/ml), and streptomycin (10 µg/ml) at 37 • C and in 5% CO 2 atmosphere (Euroclone, Italy). The culture medium is changed every 2-3 days. Cells are seeded in gridded cell-culture Petri dishes (IBIDI, 35 mm µ-Dish Grid-500). The effective area available for the cell adhesion is 40 mm 2 . The cells are stored in the incubator until the 100 % confluence condition is reached, then they are ready to be used for our experiments. Just before laser exposure, the cell medium is replaced with a DMEM solution containing a live/dead cell double staining kit consisting of calcein-acetomethoxy (AM) 0,2 µL/mL (excitation at 473 nm, emission at 490-590 nm) and ethidium homodimer-1 0,5 µL/mL (excitation at 559 nm, emission at 570-670 nm) (ThermoFisher Scientific, Italy), which are fluorescent dyes able to detect live and dead cells. In particular, the AM derivate of calcein is transported through the cellular membrane into live cells, thus colouring alive cells in green. On the other hand, the ethidium homodimer-1 is a membrane-impermeable fluorescent dye which binds to DNA, so that it enters only into dead cells, colouring them in red. We decided to insert the fluorescent probes before laser irradiation experiments and not after, as their addition requires performing some washing steps to remove the culture medium (DMEM), which often caused the detachment from the Petri dish of many cells that were damaged and stressed by the NIR pulsed illumination, thus inhibiting their count and evaluation. We note that we have not found any evidence in the literature that calcein-AM and ethidium homodimer probes may contribute to phototoxic damage. We have also experimentally verified this result, observing that, for those samples in which we managed to avoid cell detachment from the Petri dish, the damage extent in the two conditions (i.e. adding the fluorescent probes before or after irradiation) was the same.

Cell irradiation
The cells, stored in the incubator at 37 • C and in 5% CO 2 atmosphere, were quickly picked up and brought to the optics laboratory to be exposed to laser light. In the optics laboratory the temperature was kept at 27 • C to minimize thermal gradients and reduce cell stress. Irradiation experiments were performed employing an ytterbium-doped fiber laser which uses a master oscillator power amplifier (MOPA) architecture (Fidelity HP, Coherent, USA). This laser source emits pulses at 1040 nm wavelength, with repetition rate of 80 MHz, transform-limited duration 130 fs, a negative group delay dispersion variable from 0 to -120,000 fs 2 and up to 18 W of output average power. Samples were placed in a home-built inverted microscope setup equipped with white-light Köhler illumination and a CCD camera with 1024x768 pixels (Thorcam, Thorlabs, Germany), with a field of view of 200x266 µm 2 , employed to inspect the sample before the measurement and to put it in focus (by adjusting a linear translation stage along the z axis). To focus the laser beam on the sample, we employed two different optical systems: a lens with focal length f = 30 mm for the experiments presented in Section 3, using a fixed and rather large (15 µm radius, intended as half width at half maximum (HWHM)) illumination spot, and an objective (Melles Griot 10x, 0.25 NA) for the experiments presented in Section 4., based on raster scanning a smaller (6 µm radius, intended as HWHM) illumination spot. We chose these two illumination configurations with large focal spots on the sample on purpose to avoid the problem related to the vertical position of the laser beam waist with respect to the cell plane. The thin (170 µm) glass bottom of the Petri dish often presents variations in thickness. Furthermore, due to strains or bending due to gravity, it is not perfectly flat. As a result, it is rather difficult to keep the laser focus on the cells while scanning the sample over large fields of view using a high-numerical aperture objective. Using a low-NA one and a lens with long depth of focus we avoided the problem, obtaining a very good reproducibility. For the latter experiments, in particular, we combined a fast galvanometric mirror (GVS311/M, Thorlabs, Germany) with a slower motorized translation stage (8MFT-75LS05 Standa, Lithuania), moving along two orthogonal directions. In order to prevent distortion phenomena like vignetting, a 4f-relay system was implemented between the galvanometric scanner and the focusing objective.

Confocal microscopy
After laser exposure, samples were imaged using a confocal microscope (Olympus Fluoview FV10i) equipped with four diode lasers (emission wavelengths 405, 473, 559, 635 nm), an integrated incubator, and with a 10x phase-contrast objective, 0.4 NA (image size 1.2x1.2 mm 2 ). The pinhole was set to 2.5 Airy Unit. Two-colour images were acquired to detect live cells and dead cells. Images were analyzed with the open-source software ImageJ [49], counting live and damaged cells and measuring the extent of the damaged area.

Effect of fixed laser beam on cell viability
We performed a first set of experiments irradiating the sample with a fixed laser spot with 15 µm radius (see the white spot in the inset of Fig. 1) and varying both the exposure time T exp , which progressively doubled from 15 s up to a maximum exposure of 8 minutes, and the irradiance I avg , from a minimum of 0.3 MW/cm 2 to a maximum of 1 MW/cm 2 , over a 6x6 matrix of irradiated spots at 1-mm distance from each other. This translates into an average laser power between 2.26 W and 6.97 W, an energy per pulse between 28 nJ and 87 nJ and a peak power in the range of 217-670 kW. As shown in Fig. 1, the damage extends far beyond the area exposed to the laser radiation, so that the particular spatial profile of the laser spot on the sample does not play a significant role. For this reason, in analysis performed in section 3.1 and in section 3.2, we approximated the Gaussian spatial profile with a cylindrical one with the same peak value and with a diameter equal to the Gaussian FWHM. This approximation allows us to consider the irradiance as homogeneous within the laser spot. After laser irradiation, we performed a standard viability assay (as described in Section 2.1.) to determine the ability of cells to maintain a state of survival (see Fig. 1). We note that with this protocol, we are unable to detect any non-lethal damage, such as alteration in metabolism or reproductive behavior, or damage that could lead to delayed death. We see that for short exposures and/or low irradiances no damage occurred. On the other hand, we note that the damaged areas for long exposure times and high irradiances (bottom-right corner) are considerably larger than the irradiated one, and their radius increases with irradiance and exposure time.
The experiment shown in Fig. 1 was then repeated on 13 different Petri dishes, for a total of 313 irradiated sample points, corresponding to several thousands of cells, that confirmed its very good reproducibility. The average results are summarized in Fig. 2, in which we grouped the irradiance values of the laser in four bins, each of equal width (0.2 MW/cm 2 ). As already noted for the sample of Fig. 1, also for the average results no damage ever occurred for low irradiances and short exposure times (represented by the green area in Fig. 2). For higher irradiances and exposure times the damage radius increases up to an average value of hundreds of microns. Based on these observations, we first analyzed the data to determine the damage threshold as a dichotomous problem (Section 3.1). As a second step, we analyzed the subset of the damaged areas and their spatial extent, to find an empiric formula that links the radius of the circumference that encloses the dead cells to the irradiation parameters (Section 3.2).

Fig. 2.
Average radius of the damaged area for different irradiance ranges and exposure times. The green zone represents the safe zone in which no damage was detected. The red circles, whose numerical values indicate the radius, represent an indication of the damage extent.

Binary problem for estimating the damage probability
Now we characterize the damage phenomenon observed during the experiment by analyzing the data collected. Let Y be the binary variable that indicates whether or not damage has been observed on the sample irradiated by the laser (Y = 1 if cellular damage radius is greater than zero, Y = 0 otherwise). For this task, we used 301 experimental points out of 313 after outliers removal and we employed a logistic model [50], which is a sigmoidal function with monotonic increase from zero to one, commonly used in statistics as a cumulative distribution function. We computed the probability for the damage to occur as: In this formula, µ is the mean, i.e. the position of the sigmoid's midpoint for which Pr(Y = 1) = 50%, β is the logistic growth rate or steepness of the curve, which is inversely proportional to the standard deviation, and the independent variable x depends on both the exposure time T exp and the irradiance I avg with two exponents k 0 and k 1 as follows: The implemented solution aims to identify the optimal values for the exponents k 0 and k 1 , together with the mean µ and the steepness β, using the ordinary least squares (OLS) method. We iterated the fitting procedure using different combinations of k 0 and k 1 values (varying from -10 to +10 in increments of 0.1). Maximal accuracy of the prediction was found for k 0 =0.5, k 1 =1.8, resulting in µ=4.48 and β=0.98, as also reported in Table 1. The results of the fitting are shown in Fig. 3(a). The plot includes the regression line, which represents the regression equation related to the response variable written as the probability of a damage, and the experimental points. Figure 3(b) shows the same plot, but splitting the curves for different exposure times, for easier interpretation of the results.  Logistic regression can be considered as a supervised learning method belonging to the family of classification algorithms: to assign the result of the prediction to the class it belongs to (damage versus no damage), the estimated probability can be used directly as follows: where 50% is the cut-off probability. To highlight these results, in Fig. 3 we plot the experimental data points using a circle for the damaged points and a cross for those not damaged. Along this line, we coloured the background of Fig. 3 in gray for the top semi-space in which damage is probably occurring (Pr(Y = 1) >= 50%). We see that most of the circles (damaged experimental points) lie in the top gray area, while most of the crosses (non-damaged experimental points) lie in the bottom one, indicating the goodness of the model. Overall, the model presents a very high specificity of ≈95% (fraction of negative samples correctly identified as such) and a good sensitivity of ≈75% (fraction of positive samples correctly identified as such).The results are also reported in Table 1.
We can summarize these results, providing a safety rule for the microscopists employing nonlinear optical microscopy to image HeLa cells using ultrashort pulses (100∼200 fs duration) in the near-infrared (λ ≈1 µm) spectral region at high (40∼100 MHz) repetition rate. To stay on the safe side, one would like to limit the chances of damage to lower than 5%. Plugging this value into Eq. (1), we obtain x=1.475. Splitting this value using Eq. (2), we obtain the following Table 2, representing the safe working conditions. From the experimental point of view, no sample irradiated using x<=1.475 showed any damage, confirming the proper choice of the safety limit.

Continuous problem for estimating the radius of cell damage
In this second part of the analysis, we provide an empirical formula to link the radius R of the damaged area, assumed to be circular, to the exposure time T exp and the irradiance I avg . To this aim, we limited our dataset to the 155 experimental points out of 313 in which damage occurred to some extent, i.e. after removing the outliers and the green areas in Fig. 2 irradiated with low intensity and for short exposure times. Figure 4 plots the average cellular damage radius R (circles) and its standard variation (error bars) as a function of the irradiance, for various exposure times, as indicated. We observe a nonlinear dependence of R on I avg , so that we expect it to appear in the formula with an exponent larger than 1. On the contrary, we see that the exposure time must be increased by approximately one order of magnitude to double R for a given irradiance, so that we expect a rather mild dependence of the damage radius R on the exposure time T exp , with an exponent smaller than 1. We employed a regression model of this form: We iterated the fitting procedure using different combinations of the exponents h 0 and h 1 (varied from -10 to +10 with steps of 0.1). Maximal accuracy of the prediction (R 2 =0.947), using the OLS method, was found for the values of the exponents h 0 =0.3 and h 1 =1.5 (see Table 3), close to the values obtained in the previous section (k 0 =0.5, k 1 =1.8).

Heat generation and diffusion in the samples
We now investigate more in detail the type of damage itself. The observation that the damaged areas in the samples are larger than the region directly irradiated by the laser (Fig. 1), especially for long exposure times, suggests that the damage is spreading, causing the death of cells that have never been illuminated. Such behaviour can be explained in terms of the diffusion of the photoinduced heat. For this reason, in the following we consider photothermal phenomena, occurring when the main effect of the irradiation is the temperature increase. Cells undergo thermal damage when they reach temperatures around 42 • C [51]. Such high temperatures can lead to cell death or altered cell behaviour, like impaired cell reproduction [52]. In some cases, photothermal damage may be observed directly through effects such as coagulation, vaporization, carbonization and melting.
There are two main processes to consider when studying photothermal damage: heat generation and heat transport [53]. Neglecting scattering events, the former can be modelled considering that deposition of light energy in the sample follows its absorption, governed by Lambert-Beer's law. In particular, for a Gaussian beam propagating in the z-direction, the intensity within the medium is described by the following equation: with z the vertical coordinate, σ the standard deviation of the Gaussian-like profile in the radial direction (the radial coordinate being r) and α the medium absorption coefficient. To model the experimental conditions, σ has been set to be the one corresponding to a full width at half maximum of 30 µm, while absorption has been estimated by considering that in non-pigmented cells cultured in water-based media like DMEM, water constitutes the main light absorber. At our illumination wavelength (λ=1040 nm) the absorption coefficient of DMEM is α ∼ 0.2 cm −1 [54], as set in our simulations. Eventually, regarding the intensity temporal profile , I 0 (t), despite the pulsed excitation used in experiments, the high laser repetition rate (f rep = 80 MHz) allows one to simulate it as a continuous source. The value for the average intensity I 0 , which in this manuscript we call irradiance, has then been expressed as the ratio P avg /A p , with A p the beam spot area and P avg = f rep E p the average power of the train of pulses, each with energy E p . As a result of light absorption, diffusive heat transport is induced. This can be described via a heat source with volume density Q(r, z) promoting an increase of temperature across the illuminated medium. From Lambert-Beer's intensity profile, this heat source is expressed as its derivative in the z-direction, namely: according to the same notation used above, with T w the Fresnel transmission factor at the air/water interface. The additional term g(z) has instead been included in order to account for the transverse spreading of the beam due to diffraction and it acts consistently with the confinement of the beam along the direction of propagation.
Following light absorption, we study heat diffusion by solving the Fourier transport equation using a finite-element method (FEM)-based model employing a commercial software (COMSOL Multiphysics 5.6). The thermal transport problem has been implemented following standard Fourier-like formalism, namely the equation below has been solved for the temperature field T(r, z, t): Here, ρ, C p and κ are the density, heat capacity and thermal conductivity of the material, while Q is the aforementioned heat source. Values of the thermal parameters used in the simulations are referred to water after the study by Park et al. [55], and provided in Table 4. The heat transfer equation is solved across a two-dimensional axisymmetric 500 µm × 500 µm domain, by virtue of the symmetry of the sample and light excitation. To reproduce the experimental conditions, the system has been studied by defining Infinite Element (IE) domains which surround the physical domain, to simulate an infinitely extended water environment, much larger than the cell region in the sample. Beyond IEs, temperature was fixed to its equilibrium initial value, T 0 = 27 • C.
Thermal conductivity (κ) 0.6 (Wm −1 K −1 ) Absorption coefficient (α) 0.2 (cm −1 ) Figure 5(a) plots the calculated temperature, solution of Eq. (6) under the laser spot (r = 0, z = 0) as a function of time, for different irradiances. We see that within a relatively short time ( ∼500 ms, much shorter than our used exposure times) the peak temperature saturates, reaching an equilibrium value that is then maintained constant during the exposure time. As shown in Fig. 5(b), this peak temperature is linearly dependent on the irradiance and is simply the result of the dynamic equilibrium between the heat input generated by the laser absorption and its output due to diffusion. By analysing the spatial distribution of T in the radial coordinate, within less than 2 seconds the thermodynamic equilibrium is reached outside the laser spot. The obtained   temperature profile is displayed in Fig. 5(c) for three different irradiances across a region of 250 µm radius. We now employ these results to extract the values of the temperature reached by the dead cells located at the border of the damaged area (i.e. at distance from the laser spot equal to the aforementioned damage radius R). An example of this procedure is shown in Fig. 6. The first row reports the cellular damage (red cells) induced by an irradiance of 0.82 MW/cm 2 at different exposure times. We see that, at this laser fluence, for 15-s exposure there is no appreciable damage, while for longer exposure times the damage radius (highlighted in yellow in the figure) increases up to several hundreds of micrometers. The second row depicts, for the same irradiance, the results of the simulations regarding the temperature reached by the cells outside the laser spot. On each of these images, we highlighted in green an isotherm corresponding to the circumference of the damaged areas, obtained from the corresponding experimental images. The corresponding temperatures are reported below the figures. We see that for short exposure times (30 seconds) only those cells that reached temperatures close to 50 • C die, while for long exposure times (480 seconds) it is sufficient for a cell to reach temperatures around 40 • C to undergo permanent damage. We validated the results in temperature increase presented here by comparing them with those obtained with the model proposed by Macias-Romero et al. [58]. Their discrepancy is in the range of only 1-7 degrees.
This general behavior was found for all the samples and irradiances. Interestingly, we observe that the radius of the dead cells around the laser focal spot (Eq. (3)) for a given exposure time T exp only depends on the local temperature reached by the sample. To prove this, we repeated the procedure explained for Fig. 6 for all the 155 samples, retrieving the isotherm temperature at the border of the damaged areas for each irradiance. We plotted the results in Fig. 7 (blue dots), in which we grouped them as a function of the exposure time. We see that increasing the irradiance (vertical axis in Fig. 7) will cause higher temperatures on the sample, thus a larger damaged radius R, but the latter will always correspond to the same threshold temperature for cell death for a given exposure time. To highlight this behaviour, we included a red vertical line in Fig. 7 to show the mean isotherm temperature.

Methodology
While the previous set of experiments were conducted using a fixed laser illumination with a rather large spot, to simulate wide-field imaging, in this set we tested the effect of a raster scanning imaging modality on cell viability. Some studies [59,60] already proposed various power and irradiance thresholds to avoid tissue damage, and our goal here is to build on those studies to cover a wider range of parameters.
To this aim, we employed the combination of galvanometric mirror and motorized translational stage on two perpendicular axes, as described in Section 2.2. The laser beam followed a serpentine pattern to cover the whole sample area, changing direction after each horizontal scan is completed. In particular, we systematically varied the irradiance between 0.9 MW/cm 2 and 6.2 MW/cm 2 , the scanning speed (0.25 mm/s, 0.5 mm/s and 1 mm/s) and the number of scanned lines (20, 50 and 200) within the same illumination area of 50x50 µm 2 . In this way, the vertical distance between two consecutive horizontal scans is inversely proportional to the number of scanned lines. Given the irradiance, the laser repetition rate and the numerical aperture of the objective described in section 2.2, this translates into an average laser power between 1 W and 7 W, an energy per pulse between 12.5 nJ and 87.5 nJ and a peak power ranging from 96 kW to 673 kW and a peak power density between 85 and 595 GW/cm 2 . Considering the scanning beam diameter (d), 12 µm, the number of lines (N) and the vertical dimension of the illumination area, 50 µm, it is possible to define the vertical line density (VLD) as: With this parameter we consider the fact that scanning lines overlay during scanning. The energy density, also called light dose, is defined as irradiance (I avg ) * exposure time (T exp ). As exposure time, in this experiment we considered the total time the single cell is illuminated by laser light, and we calculated it in this way: where A spot is the area of the scanning laser spot and v s is the scanning speed. As a result, the energy density ranges from 0.04 to 11.2 MJ/cm 2 .

Cell counting
After laser irradiation, the samples were placed in the confocal microscope with incubator described in Sec. 2.3, where images of the scanned areas were taken just after the laser irradiation and after 18 hours. This enabled us to examine the onset of cell death over a longer period and obtain a more reliable result, as radiation-induced apoptosis in some cells may not be immediately visible. In Fig. 8 we report the image of one of the samples (18 hours after the irradiation), showing 54 different irradiated areas varying the irradiance (vertical axis), the scanning speed and the number of scanned lines (horizontal axis). Each 50x50 µm 2 area contained an average of 10 cells. We repeated this experiment on three separate samples, obtaining consistent results across more than 1600 cells. We only counted cells in the area actually irradiated by the laser, even if sometimes the damage exceeded this area, such as in the cases of 200 scanned lines lines (third, sixth and ninth columns in Fig. 8) irradiated with the maximum irradiance (last row in Fig. 8).

Data analysis
The goal of this analysis is to find a data-driven law that traces the relationship between the percentage of damaged cells and the experimental variables, i.e. laser power (P avg ), scanning speed (V s ) and number of scanned lines (L). Let Y be the percentage of damaged cells over the total in the scanned area. Similarly to the previous section, we modeled Y as a logistic function Y = f(x) (see Eq. (1)), where the independent variable x is dependent on the experimental variables with given exponents k 0 , k 1 and k 2 : In particular, we identified the optimal values for the exponents k 0 , k 1 , k 2 for the entire dataset using 187 samples (no outliers detected) by iteratively fitting the data using different combinations of values in the -10 to +10 range (in increments of 0.01). Using the OLS method, the best results that minimize the sum of the squares of the differences between the observed data and those of the curve representing the function itself are reported in Table 5. Figure 9(a) plots in black the fitted logistic function with data collected 18 hours after the laser irradiation. We also added in red the result using data collected just after the irradiation. The black curve is left-shifted, meaning that the number of damaged cells increased, and this confirms the inaccuracy of the data collected immediately after irradiation. Only using data collected after 18 hours, to better appreciate the goodness of the procedure for the different number of scanned lines, we also show in Fig. 9(b) a separate scatter plot for 20 (green), 50 (yellow) and 200 (cyan) scanned lines, with the same logistic function, in which the abscissa is now just a function of the average power and the scanning speed. In this way, we clearly see that increasing the number of scanned lines in the same field of view induces higher damage on the cells because they are illuminated several times during the imaging process. Fig. 9. a) Logistic functions F(x, µ, β) with x = P k 0 avg · L k 1 · V k 2 s , in red obtained with data collected just after the irradiation and in black with data collected after 18 hours. b) Logistic functions F(x, µ, β) with x L = P k 0 avg · V k 2 s for each value of scanned lines L, with data collected after 18 hours.
Looking at the optimal values reported in Table 5, we see that the dependence on P avg is almost linear (k 0 = 1.0), while there is an inverse dependence on the scanning speed (k 2 = -0.08): as expected, the percentage of damaged cells decreases as the speed increases. However, the exponent is very close to zero, as also confirmed looking at the heatmaps in Fig. 10: the percentage of damaged cells changes much more noticeably when the laser power varies (vertical axis) than when the scanning speed varies (horizontal axis).
As done in Section 3 for fixed laser illumination, valid for wide-field imaging configuration, we provide here an empirical safety rule for the parameters in the raster scanning illumination configuration to keep the damage probability below 5%. In this case, the safety threshold corresponds to x = 5.367. Only one damaged experimental point was reported below this threshold, thus confirming the goodness of the 5% safety limit. Table 6 reports the safe working conditions, splitting the x values using Eq. (9).

Discussion and conclusions
In this work, we explored the photodamage mechanisms induced in live HeLa cells by NIR ultrashort laser pulses, conducting two experiments with different laser exposure methods. In the first experiment, the laser was held fixed on the sample, while the second experiment was conducted under raster-scanning configuration. There is a well known relationship between photodamage rate R and laser power applied: R ∝ P n p f τ, where P p is the laser peak power, n is the order of photodamage mechanism, f is the laser repetition rate and τ is the pulse duration [29,61,62]. To analyze our data, since we decided not to vary repetition rate and pulse duration, we chose to evaluate the dependence of damage on average power, without losing generality. In the first experiment we used the irradiance, I avg = P avg / spot area [MW/cm 2 ], which is equivalent to average power since we kept the spot size fixed for the duration of the experiment.
Considering the result of the 'Binary problem for estimating the damage probability', Subsection 3.1, we obtained 1.8 as exponent for irradiance, while in the 'Continuous problem for estimating the radius of cell damage', Subsection 3.2, we obtained an exponent of 1.5. In both cases, the exponent is < 2, and this means that the photodamage is mainly caused by a first-order mechanism, therefore photothermal, with some minor contribution from higher-order processes. This conclusion is confirmed also in Subsection 3.3, in which, using a thermodynamic model, it is demonstrated that potentially lethal temperatures are reached in correspondence with the damaged areas. Regarding the contribution from higher-order processes, photochemical damage with ROS formation may be present, also enhanced by heat shock [63,64]. It must be noted that another mechanism comes into play when cells are stressed: the exchange of signals between them. Apoptosis can either be started by the cell itself (intrinsic pathway) or it can be triggered by signals from other cells, that alert their neighbours of a "dangerous" situation (extrinsic pathway) [65]. This might explain why some cells that are not exposed to particularly high temperatures (40.0 • C, see bottom-right panel in Fig. 7) still start the apoptosis process, mainly for long irradiation exposure (480 s in this case): the infrared laser radiation generates heat, that acts as environmental stressor. In addition, it has been demonstrated that heat shock triggers the death receptor apoptotic pathway and, as already said, ROS production, leading to cell death [66][67][68]. In addition, considering the theory of Brownian motion, random walk and diffusion, we know that diffusing molecules travel over a distance whose square value increases on average linearly with time [69]. Thus, this hypothesis of exchange and diffusion of signals between cells is also confirmed by the fact that we observed that the radius of the damaged cells depends on the exposure time raised to a power close to 0.5 (see section 3.2).
The main finding of this set of experiments is therefore the fact that cells can withstand high temperatures for a short amount of time, while they will die if exposed for longer periods of time to even mild temperatures. This is a valuable piece of information for microscopists performing nonlinear imaging with ultrashort laser pulses. It provides a clear proof that it is in general more convenient to illuminate the samples with high irradiances: indeed, they provide strong signals (which scale like the second or third power of the irradiance for nonlinear processes such as second-harmonic generation, two-photon excited fluorescence or third-harmonic generation) so that the sample can be imaged in a short amount of time without harm. On the other hand, employing low irradiances would call for much longer imaging times, thus risking to create permanent damage to living cells.
Moving on to the second experiment, the one with the scanning laser beam, we obtained an exponent equal to 1.0 for the average laser power. In this case the linear damage component is even more predominant, and this can be explained considering the different exposure times applied. In the first experiment, the exposure time ranged from 15 s to 480 s, while in the second experiment it ranged from 45 ms to 1.8 s, too short to induce the formation of ROS in this parameters range.
Furthermore, the data-driven approach used to analyze the data, provides a valid tool to evaluate the experimental conditions to be used to avoid damaging the sample, if the same range of parameters are used. The 5% damage probability appears to be a proper safety threshold parameter. In fact in our experiments, photodamage was found only once below this threshold.