Optimizing ultrafast illumination for multiphoton-excited fluorescence imaging

We study the optimal conditions for high throughput two-photon excited fluorescence (2PEF) and three-photon excited fluorescence (3PEF) imaging using femtosecond lasers. We derive relations that allow maximization of the rate of imaging depending on the average power, pulse repetition rate, and noise characteristics of the laser, as well as on the size and structure of the sample. We perform our analysis using ~100 MHz, ~1 MHz and 1 kHz pulse rates and using both a tightly-focused illumination beam with diffraction-limited image resolution, as well loosely focused illumination with a relatively low image resolution, where the latter utilizes separate illumination and fluorescence detection beam paths. Our theoretical estimates agree with the experiments, which makes our approach especially useful for optimizing high throughput imaging of large samples with a field-of-view up to 10x10 cm2.


Introduction
Imaging utilizing two-photon excitation of fluorescence (2PEF) and, since recently, threephoton excited fluorescence (3PEF), has significantly enhanced biological microscopy by improving image quality and acquisition speed, reducing photo-damage to the sample, as well as by allowing deeper sample depth penetration [1][2][3][4][5][6][7][8][9][10][11]. While a broad variety of new imaging modalities continue to be demonstrated that take advantage of the high peak pulse intensity of different femtosecond laser sources, optimizing the illumination and detection conditions still poses an issue. One reason for this is because the average rate of multiphoton excitation depends not only on the average incident power, as in the case of linear (one-photon) excitation, but also on the instantaneous flux of photons, i.e. peak intensity of the pulses. This imposes constraints on the photon budget, which, along with such factors as the sample damage threshold, maximum laser power available, and various sources of experimental noise, makes finding the optimal illumination-and detection conditions a challenging task.
In this paper we provide quantitative analysis of the maximum imaging rate that can be achieved under realistic conditions in femtosecond multiphoton excited fluorescence systems. Faster imaging is imperative for capturing real-time processes such as the activity of individual nerve cells, which occur on a sub-millisecond time scale, and is also critical for reducing long scan times that currently limit multi-photon imaging of large-area samples. We achieve our goal by taking into account different illumination and focusing conditions, sample properties, and the signal-to-noise ratio (SNR) of the image data, where the latter may be affected by photon shot noise, laser pulse energy fluctuations and fluorescence detection noise. As the first step, we consider a generic multiphoton microscope arrangement, with coaxial illumination-and fluorescence detection beam paths, using the tightly-focused beam of a 100 MHz pulse rate femtosecond oscillator. We evaluate the maximum achievable shot noise limited rate of imaging depending on the two-photon absorption (2PA) and three-photon absorption (3PA) cross section values of the fluorophores and other key parameters such as average laser power, sample damage threshold etc. Some applications such as early cancer detection [12] and high throughput screening [13] require much larger samples, that cannot be addressed by conventional coaxial illumination setups. For this purpose, we model a wide field-of-view (FOV) imaging system, where the illumination-and the fluorescence detection beam paths are separated, and where the illumination uses a ~kHz pulse rate amplified femtosecond laser. In order to verify our model, we compare the theoretical results with published data about different multiphoton microscope systems. We also perform experiments, where we determine the maximum imaging rate in a CCD camera based, ~10x10 cm 2 FOV, 2PEF imaging system using different 1 kHz amplified laser sources, and show that our theoretical analysis has a good correlation with the experiments. To our best knowledge, this paper presents the first quantitative analysis of femtosecond illumination for wide FOV 2PEF imaging as well as the first theoretical treatment of imaging speed optimization for multi-photon fluorescence microscopes. Figure 1(a) shows the femtosecond laser beam that emanates from the output aperture of a high-NA microscope objective and converges onto a diffraction-limited 3D spot inside a sample of thickness, d. If we would neglect potential side-effects associated with the femtosecond illumination, such as photo-bleaching and saturation of the 2PA transition, then the number of the fluorophores excited by a single laser pulse may be expressed as Eq. (1) [3,9,14]

2PEF and 3PEF imaging with coaxial illumination
where σ 2PA is the 2PA cross section (in units of GM = 10 −50 cm 4 s photon −1 ) of the fluorophores at the illumination wavelength, λ L , I L is the photon flux (number of photons per unit area and per unit time), C is the concentration of 2PEF-active fluorophores, and the integration is carried out over the duration of the pulse and over the entire volume of the sample. If we assume, for simplicity, that the fluorophore concentration is constant throughout the sample, then most of the 2PEF signal will be produced by the volume confined inside the focal region (voxel), where the photon flux is the highest. In the case of a Gaussian spatial beam profile, the voxel may be defined by the closed surface, where the peak photon flux decreases by factor e −2 from its maximum value, where ω 0 is the beam waist related to the numerical aperture of the objective, ω 0 = λ L (π NA) −1 , and the length is given in meters. If we approximate the temporal intensity profile with a Gaussian with FWHM value, τ p , and if the Rayleigh length of the focus, Δz R = (πω 0 2 )/λ L , is much less than the sample thickness, Δz R << d, e.g. as shown in Fig. 1(a) where P ave is the average excitation power in W, g is the pulse repetition rate in s −1 , h is Planck's constant and c is the speed of light (both in MKS units). The above relation illustrates the well-known fact that the number of excited fluorophores is proportional to the square of the incident power and to the 2PA cross section. It also shows, provided that the voxel is entirely contained inside the sample, the independence of the excitation rate on NA, or on how tightly the beam is focused. Equation (3) may be recast in a more familiar way, that was first introduced by Denk et al. [3]: where the quantity enclosed in the square brackets is the number of photons absorbed per fluorophore per pulse. For multiphoton microscopy the focus point is scanned in the x-y directions until it covers the entire FOV area, e.g. in a raster pattern shown in Fig. 1(b), and a 2D image is built pixelby-pixel from the recorded 2PEF signals. The lateral resolution is determined either by the focus spot size or by the raster step size, whichever of the two is larger. Here, for simplicity, the 2D focus spot is shown as a square pixel of area, ΔS = Δx 2 , and the raster-scan step size is set equal to Δx. The average number of the fluorophores excited in the n th pixel is then proportional to the number of laser pulses delivered at that pixel: where τ dw is the pixel dwell time. The total number of pixels in the FOV is, The inverse of the minimum time required to scan the whole FOV is called the frames-persecond, or FPS, and may be expressed as: 1 . In the coaxial arrangement, the objective used for illumination also collects the fluorescence photons, which are then passed on to a photo-detector (not shown in Fig. 1). The number of fluorescence photons detected from one laser pulse is: where ε eff is the fluorescence detection efficiency and η MO is the fluorescence photon collection efficiency. Assuming that the fluorescence photons are not scattered by the sample, and using the definition of the numerical aperture; NA = nsin(θ), where n is the index of refraction of the sample and θ is the half angle of the cone of light collected by the objective, we can express the fluorescence collection efficiency as [15]: Spectral band-pass filters are used to reject scattered laser light as well as any other spurious photons, while allowing the fluorescence photons to pass through to the detector. The overall fluorescence detection efficiency may be expressed as: where Φ FL (λ F ) is the differential fluorescence emission quantum efficiency of the fluorophores, φ F (λ F ) is the spectral filter transmission function and κ D (λ F ) is the quantum efficiency of the photo-detector. The integration in Eq. (10) is performed over all fluorescence wavelengths. Ultimately, the utility of any image depends not only on the highest attainable resolution and fastest FPS, but also on whether the sought-after information may be extracted in the presence of the inevitable noise. In order to quantify the latter, we introduce the per-pixel SNR defined as: where δ shot stands for the photon shot noise, δ laser is the noise due to laser pulse-to-pulse energy fluctuations and δ el represents detection noise, e.g. due to electronic amplifiers. The pulse-to-pulse energy fluctuations of mode-locked femtosecond lasers operating at high pulse rates, g > 1 MHz, are usually small, while δ el can be minimized by using photon counting detection schemes. This leaves δ shot as the main contributing factor to the noise, in which case, Eq. (11) may be expressed as: By solving Eq. (5) for τ dw , and using the relations Eqs. (8), (12), and (3), we can express Eq. (7) for the maximum attainable frame rate as: Equation (13) illustrates the well-known fact that the highest attainable imaging rate in a 2PEF microscope system is limited by the maximum incident power, i.e. excitation photon budget [7]. Figure 2(a) presents the value of FPS max , plotted as a function of P ave , for g = 100 MHz, SNR min = 10, λ L = 790 nm, and FOV = 100x100 μm 2 . Other parameter values used in this calculation are listed in the figure caption. Increasing the average power allows for a decreased dwell time, which in turn increases the frame rate. However, damage to the samples usually limits the average laser power to about, P ave < 50 mW. On the other hand, if the fluorophores possess large σ 2PA , then the 2PA transition may approach saturation. To obtain a quantitative estimate of the saturation effect, we numerically solve the kinetic equations for the excited state population in the two-level approximation [16,17]: where ρ f is the probability of finding the chromophore in the excited state and γ f0 is the spontaneous excited state population relaxation time back to the ground state, 0. Because in our case the excitation pulse duration is much less than γ fo −1 , Eq. (14) predicts that a maximum of half of the fluorophores can be excited during one laser pulse, which prompts us to set the saturation limit at half of the maximum value. In Fig. 2(a) the two-photon saturation limit is represented by the horizontal dashed line. With all these factors taken into account, we obtain, FPS max ~10 3 . At this point it is useful to introduce the rate of imaging, right vertical axis of Fig. 2(a), defined as, Im .
The estimated maximum attainable rate is then, ~10 −1 cm 2 s −1 . Instantaneous three-photon absorption (3PA) is a process where the transition from the ground to the final excited state occurs due to simultaneous absorption of three longerwavelength photons. The version of multiphoton microscopy using three-photon excitation of fluorescence (3PEF) takes advantage of the ability of longer-wavelength photons to penetrate even deeper into scattering tissues [23], and can also potentially achieve a tighter excitation confinement than 2PEF. In the case of degenerate 3PA, the 3PA equivalent of Eq. (1) is: where σ 3PA is the 3PA cross section at the wavelength λ L . Integrating over the voxel volume gives the number of fluorophores excited by one laser pulse: If we assume, as before, that the accuracy of measuring the 3PEF signal is limited by the photon shot noise, then the corresponding maximum attainable frame rate is: Peak values for the three-photon cross section of organic fluorophores are rather small, typically in the range, σ 3PA ~10 −82 −10 −84 cm 6 s 2 photon −2 [11,24,25]. In practical terms, this means that the relative efficiency of 3PEF is lower compared to that of 2PEF, assuming similar illumination conditions. On the other hand, the maximum average power is still limited by the sample damage threshold, P ave < 50 mW. A suitable solution to this issue was demonstrated by Xu et al. [11], who used a femtosecond laser with about two orders of magnitude lower pulse repetition rate, which boosted the peak photon flux by the same factor. Figure 2(b) presents the dependence of FPS max on P ave for g = 1 MHz, and for different σ 3PA values. The saturation limit is represented again by the horizontal dashed line. Other parameter values are listed in the figure caption. The estimated maximum imaging rate, ~0.6x10 −2 cm 2 s −1 , turns out to be about one order of magnitude less compared to the 2PEF rate at the same average power. Unfortunately, there is a lack of commercially available, near-IR few-MHz pulse rate femtosecond lasers providing sufficient output power. The symbols in Fig. 2

2PEF and 3PEF imaging using separated illumination-and fluorescence beam paths
Our goal here is to perform a similar photon budget and maximum imaging rate analysis with a much larger FOV. For example, if the sample encompasses a 10 cm diameter Petri dish, then a straightforward extrapolation of the above estimates predicts that it would take at least 10 3 s to acquire a two-photon image of the whole plate. Even if we would disregard the excessively long acquisition time, aberrations caused by large beam scan angles [26,27] would render the coaxial illumination arrangement rather incompatible with large objects. An alternative approach is to increase N 2PA by increasing the peak photon flux. This may be achieved using amplified femtosecond lasers that operate at ~kHz pulse rates. Such lasers have been commercially available since the 1990-ies and are capable of delivering up to 10 6 times higher peak power compared to mode-locked femtosecond oscillators. However, to take full advantage of this high power, an entirely different illumination-and detection system design strategy is on order. First of all, diffraction-limited focusing becomes impractical because the resulting extremely high peak intensity may damage not only the sample, but also the objective along with other setup components. Secondly, with a loosely focused illumination beam path, it becomes advantageous to separate it from the fluorescence collection beam path because the former no longer needs high NA optics. Figure 3(a) shows the principle of a wide FOV imaging setup, where the multiphoton excitation is induced by a loosely-focused amplified femtosecond laser beam, while the fluorescence is collected by a camera objective. Due to the small effective NA of the illuminating beam, the corresponding Rayleigh length is now comparable to or even exceeds the sample thickness, Δz R > d. As will be discussed below, the shape of the beam will be optimized according to the available laser power and properties of the sample, but its lateral size rarely exceeds ~0.1 -2 cm. To cover the entire FOV, the beam is moved in the x-y directions, e.g. in a raster scan pattern depicted in Fig. 3(b), where we have again assumed, for simplicity, a rectangular beam area, ΔS = ΔxΔy, and that the scan proceeds in equidistant steps without overlap between the neighboring positions. The fluorescence is collected by a lens, as shown in Fig. 3(b), and then focused on the imaging sensor of a camera (not shown). The image resolution is now determined by the lens and the camera sensor, while the maximum FOV is determined either by the camera FOV or by the maximum beam scan area, whichever of the two is smaller. For simplicity, let us assume that the resolution equals the pixel size, {Δx pix , Δy pix }. where ε eff is the camera detection efficiency and η CO is the fluorescence collection efficiency. When considering the image noise we need to account for the fact that the output of femtosecond amplifiers is noisier than the output of mode-locked oscillators. For example, a typical regenerative Ti:Sapphire amplifier displays short-term root mean square pulse-to-pulse energy fluctuations of about 5 -20% (RMS laser = 0.05 -0.2). Let us assume, for simplicity, that the photon shot noise and the electronic camera detection noise are small and may be neglected compared to the random laser pulse energy fluctuations. In this case, the image pixel SNR increases in proportion to the square root of the number of laser pulses that are averaged at each pixel, which leads us directly to the following simplified expression for the maximum imaging rate:  Fig. 2(a). Because these practical ROIm values fall substantially short of the estimation, we take this as a clue that a rectangular beam shape may not be accurate enough. The excitation is better described by Gaussian functions, which in our case may be expressed as: where Δx and Δy are the beam spot sizes. It is also important to consider that the illumination beam position is moving continuously over the FOV, rather than in discrete, non-overlapping steps. In practical terms, this means that each pixel is being illuminated not just by one pulse but by multiple pulses, with a different effective photon flux every time, as the spatial distribution, given by Eq. (22), moves over the corresponding area. Furthermore, because the number of photons available per illumination pulse is limited by the maximum power of the laser, it is imperative to optimize the beam size (ΔS = ΔxΔy) as well as the scan speed (dwell time), such that, for each image pixel, the desired minimum SNR is achieved. In other words, a realistic simulation requires an illumination scheme with a moving, low pulse repetition rate Gaussian beam, where each pulse is added up and averaged individually. These rather involved, inter-dependent conditions may be summarized in the following relation: ., where the empirical parameter β s represents the optimal average power density that is to be determined for each particular experiment and/or type of sample. All the above aspects are taken into account by our numerical model that is described in detail in the Appendix. The results of the numeric simulation are presented in Fig. 4(a). The calculation assumes the following parameter values: σ 2PA = 10 GM, g = 1 kHz, τ p = 150 fs, C = 2x10 21 m −3 , η CO* ε eff = 1x10 −4 , and SNR min = 10. The power parameter is, β s ~2.0x10 14 W m −2 , which corresponds to the optimal illumination intensity for imaging of E. coli bacteria colonies expressing fluorescent protein mutants [13]. As before, the lower horizontal axis of the figure corresponds to the average power, which is in this case limited to 4W. The upper horizontal axis shows the corresponding optimal beam spot area determined from Eq. (23). As expected, the maximum imaging rate depends strongly on the laser noise and on the available maximum power, and in the case of RMS laser = 0.2, with σ 2PA ~10 GM peak cross section fluorophores, reaches, ~0.01 -0.1 cm 2 s −1 . According to this calculation, the minimum number of pulses needed to achieve SNR = 10, is 50, 200, or 800 pulses for RMS laser ~5%, 10%, and 20%, respectively. We may conclude that given a typical amplified femtosecond laser, a 2PEF image may be collected from an entire 10 cm diameter Petri dish within less than an hour. Typical scientific-grade cameras have a read-out noise of about ~10 −100 counts per pixel. If we would assume that the laser noise is absent, then, to guarantee a minimum per-pixel SNR = 10, the fluorescence signal should be at least ~10 2 -10 3 counts per pixel. Note that at such a relatively high signal level the shot noise may be disregarded. However, a notable complication arises if the distribution of the fluorophores in the sample is not uniform. Figure  5(a) shows a representative case, where the fluorophores are concentrated primarily in distinct features, such as bacterial colonies expressing fluorescent proteins [13]. The goal of the multiphoton imaging experiments is, in this case, to determine the total fluorescence emitted by each feature, which are randomly distributed over the entire FOV. Therefore, we need to augment our model by including the illumination of each feature. Figure 5(b) provides the x-y view of the model sample, along with the outline of the asymmetric Gaussian beam and the corresponding raster scan pattern. The fluorescent features are depicted as circles, where the radius is larger than the image pixel size, but less than the illumination spot size. If the sample is homogeneous, or if the pulse repetition rate is high, then the ROIm would not need to account for the details of the illumination pattern. In our case, however, we need to calculate explicitly how much fluorescence is emitted by each feature excited with each pulse. To facilitate this, we divide the 10x10 cm 2 FOV into an illumination matrix comprised of i max elements in the x direction and j max elements in the y direction, where each element {i, j} represents one excitation beam position. The number of fluorophores excited in the k th feature may then be then expressed as: where the photon flux I L (x, y) is given by Eq. (22) and C k (x,y) is the effective fluorophore concentration in the k th feature. Details of the numerical simulation are given in the Appendix. Briefly, we start by creating a random pattern of fluorescent features to be detected, and determine the optimal beam area, ΔS, corresponding to the average power value, P ave , using Eq. (23). We set the Gaussian beam parameters such that ΔxΔy = ΔS, and the degree of asymmetry matches the experimental conditions. Because the FOV is constant, the choice of i max and j max determines the raster step size in the x-and y-direction, respectively. We now calculate the total number of fluorescence photons detected from each feature using Eq. (24) and (20) for a range of different scan speeds, or equivalently, different imaging matrix sizes, which allows us to determine when the minimum SNR per feature is achieved. The results are shown as the blue line in Fig. 4(a), and the corresponding parameter values are listed in the figure caption. At this point we may conclude that, given an ideal, noise-free 1 kHz amplified femtosecond laser, one should achieve rather fast 2PEF imaging of large area samples at a rate of at least ROIm ~1 cm 2 s −1 , even at moderate average powers < 0.1 -1 W.
It is also of interest to estimate how a 100 MHz source would perform under similar conditions. This simulation is shown as the dashed line in Fig. 4(a). Even though relatively high average power > 1W will be needed to achieve ROIm ~0.1 cm 2 s −1 , this type of illumination source may also be potentially considered, provided that issues associated with the requirement of rapid scanning of a tightly-focused beam over a large area without loss of the beam position accuracy can be resolved. Figure 4(b) shows the maximum imaging rate obtained by applying the same calculation to three-photon excitation of fluorescence with different σ 3PA values. Since a further increase of the peak photon flux per pulse is not practical, useful imaging rates may be achieved only for fluorophores with exceptionally large 3PA cross sections, on the order of 10 −80 cm 6 s 2 or more. Figure 6 shows the principle of the 2PEF imaging experimental setup, which was described, in part, in [13,28]. The output beam of a 1 kHz pulse rate femtosecond laser is shaped by two lenses (L1, L2) and is scanned over the area of a standard 10-cm diameter Petri dish using a bi-axial motorized mirror mount (Zaber T-OMG). The fluorescence emitted from the sample is collected with an objective lens and is detected with a TE-cooled CCD camera (Hamamatsu C4742-98) with an estimated image resolution of 200 μm. A stack of band-pass filters (F1) placed in front of the camera objective cut off scattered laser light. A reference photodiode (PD) was used to monitor and correct for long-term (~minutes) changes of the laser power. As the model sample, we use E. coli colonies expressing EGFP-type fluorescent proteins [13,28]. The beam steering and data collection were controlled by a LabView routine, and the data was analyzed using custom MatLab code as described in [13]. The symbols in Fig. 4(a) show the maximum imaging rate obtained with three different 1 kHz femtosecond lasers sources: (a) Clark MXR CPA-1000 Ti:Sapphire regenerative amplifier (solid circle); (b) Coherent Legend-HE Ti:Sapphire regenerative amplifier (empty rectangle) and (c) Light Conversion TOPAS-C optical parametric amplifier (empty circle). In both regenerative amplifier laser systems, the pulse-to-pulse instability was, RMS laser ~0.1. Because the maximum average power of the first system was limited to P ave ~500 mW at λ L = 790 nm, the resulting maximum imaging rate was about ~0.08 cm 2 s −1 , while the maximum average power of the second laser system, P ave ~1300 mW, allowed, at the same excitation wavelength λ L = 790 nm, about twice as fast imaging rate, ~0.2 cm 2 s −1 . In the third case, the wavelength was converted by OPA to λ L = 925 nm, which is closer to the peak of the 2PA spectrum of EGFP than λ L = 790 nm [13]. However, because the OPA had less maximum output power, P ave ~80 mW, and also because this source had more pulse-to-pulse instability, RMS laser > 0.15, the maximum imaging rate was lower, ~0.007 cm 2 s −1 . Figure 7 shows plots of the corresponding beam profiles. The beam area was determined in all three cases according to Eq. (23), but the aspect ratio, i.e. Δx vs. Δy was chosen based on convenience of steering the beam over the sample.

Conclusions
We have analyzed the performance of femtosecond multiphoton excited fluorescence imaging systems using different illumination-and focusing modalities, while taking into account the minimum required SNR and dominant sources of noise, such as photon shot noise, laser pulse energy fluctuations and fluorescence detection noise. In the case of a tightly-focused beam in the coaxial arrangement, using a 100 MHz illumination source to produce 2PEF, the maximum rate of imaging is limited by the fluorescence photon shot noise, by light-induced damage to the sample, and by saturation of the two-photon transition, and may reach about ~10 −1 cm 2 s −1 , depending on the 2PA cross section and concentration of the fluorophores. In the case of coaxial 3PEF imaging with a 1 MHz illumination source, the effect of saturation limits the maximum imaging rate to about ~10 −2 cm 2 s −1 . However, due to aberrations, in a coaxial illumination system it is difficult to access a FOV much larger than a few mm 2 . To achieve wide FOV 2PEF imaging, up to 10x10 cm 2 , one can take advantage of ~kHz pulse rate femtosecond amplifiers and use separated illumination and fluorescence imaging beam paths. We show that the highest rate of imaging, up to 1 cm 2 s −1 , may be achieved with a virtually noise-free amplified laser, and is limited just by the camera detection noise. If the laser pulse-to-pulse fluctuations are significant, then the maximum imaging rate decreases due to time averaging of the signal, leading to a maximum ROIm ~0.1 -1 cm 2 s −1 . We perform experiments using three different 1 kHz pulse rate femtosecond lasers with different characteristics and find good agreement with our theoretical estimates.
Our results may be used to optimize multiphoton imaging systems that are designed to perform challenging tasks such as high throughput screening of two-photon fluorescent proteins, detection of early-stage cancer, or real-time monitoring of the brain cell activity in behaving animals. Similar analysis to that performed here could benefit other novel imaging modalities such as photoacoustic imaging [29].