Focus image scanning microscopy for sharp and gentle super-resolved microscopy

To date, the feasibility of super-resolution microscopy for imaging live and thick samples is still limited. Stimulated emission depletion (STED) microscopy requires high-intensity illumination to achieve sub-diffraction resolution, potentially introducing photodamage to live specimens. Moreover, the out-of-focus background may degrade the signal stemming from the focal plane. Here, we propose a new method to mitigate these limitations without drawbacks. First, we enhance a STED microscope with a detector array, enabling image scanning microscopy (ISM). Therefore, we implement STED-ISM, a method that exploits the working principle of ISM to reduce the depletion intensity and achieve a target resolution. Later, we develop Focus-ISM, a strategy to improve the optical sectioning and remove the background of any ISM-based imaging technique, with or without a STED beam. The proposed approach requires minimal architectural changes to a conventional microscope but provides substantial advantages for live and thick sample imaging.


Working principle of Focus-ISM
The point-spread function (PSF) of the scanned-image generated by the element at position x d of a detector array is given by the PSF of a laser scanning microscope, just with a shifted pinhole aperture [8,9]. The equation is where x s = (x s , y s ) are the scan coordinates on the sample plane in focus, p is the pinhole aperture that describes the sensitive area of the element of the detector array located at x d = (x d , y d ), and h exc and h det are, respectively, the excitation and detection normalized PSFs ( h exc = h det = 1). For STED-ISM, we also have to consider the reduction of the effective fluorescence region, which is equivalent to a reduction of the effective excitation PSF where h (0) exc is the confocal excitation PSF, and η is the depletion function, which describes the reduction of fluorescence emission under exposure of the depletion beam. ς is the saturation factor, defined as (S3) ς = I STED /I sat where I STED (x s ) is the doughnut-shaped intensity distribution of the STED beam and I sat is the saturation intensity, namely the intensity of stimulated emission needed to have spontaneous emission as likely as stimulated emission. In other words, the spontaneous emission rate k F equals the stimulated emission rate k SE when I STED = I sat . In order to calculate the saturation factor ς using a pulsed laser, we need to take into account the dynamics of the fluorescence intensity F (t) under stimulated emission [8]. The fluorescence decay is where T is the pulse duration of the STED beam. Assuming a time-gated detection in the microscope (namely, the light is collected only after the STED beam pulse is over), the depletion function for the gated pulsed STED implementation reads Finally, the equation for the depleted excitation PSF is Since ς(x s ) is doughnut-shaped, the excitation PSF is more depleted at its borders and its effective width is smaller. Notably, the STED excitation PSF is no longer normalized.
For the rest of the manuscript, we will indicate with h exc a generic excitation PSF with effective width depending on the saturation factor. Thus, the scannedimage generated by the detector element at position x d can be written as where o(x s ) is the specimen, here considered as a distribution of light emitters.
In this work, we consider the limit case where each element of the SPAD array can be described as a shifted Dirac delta function. Substituting p(x s ) with δ(x s − x d ) in equation S1 we obtain In order to grasp the core idea behind image scanning microscopy (ISM), we first introduce its concept using a Gaussian approximation and later proceed with the general case.

Gaussian model
In this subsection, we simplify the model approximating the emission and excitation PSF with a Gaussian distribution with circular symmetry, where µ is the mean and σ is the standard deviation. Thus, the resulting PSF is Notably, the product of two Gaussian function is another Gaussian function with the following properties: where µ(x d ) is known as shift-vector, and σ ism is independent of x d . The intensity of the resulting Gaussian function is also rescaled by the following rescaling factor Thus, the resulting PSF can be written as Neglecting the Stokes-shift, we can assume σ exc ∼ σ det to obtain the well-known √ 2 gain in resolution. However, the images generated by each detector element are shifted with respect to the central element. Indeed, the scanned-image is If all the images are simply summed together, the resolution gain is lost, and the resulting image is equivalent to that generated by a traditional confocal microscope with an open pinhole (as large as the detector size). The pixel reassignment (PR) algorithm shifts back each image and constructs a new image summing the reassigned images thus conserving the gain in resolution and enhancing the signal-to-noise ratio (SNR).

General model
The typical assumption behind ISM is that the complete PSF can still be written as a new shifted function Thus, the scanned-images generated by each detector are all of similar shape, but shifted, Instead of using the theoretical values found using the Gaussian approximation, the shift vectors are found by the adaptive pixel reassignment (APR) algorithm. This latter first calculates the correlogram with the purpose of making all the scanned-images identical, except for an intensity scale factor Rather than considering the scanned-image, the ISM image formation process can be seen equivalently from the perspective of the detector. In fact, the detector array can be considered as a small camera, capable of acquiring widefield images that we call micro-images. The equation of these latter can be found as follows: As expected, the micro-image is simply given by the wide-field image formation law, but applied to the illuminated object, i.e. the object multiplied by the excitation PSF. After APR, the reassigned micro-image is For clarity, we now rewrite the new coordinates of the micro-images as where we have used the theoretical value of the shift-vectors found using the Gaussian approximation for simplicity. In this new coordinate system, the post-APR micro-image is In other words, the effect of APR is a digital shrinking of the excitation and detection PSFs with respect to the object. In detail, the standard deviations of the PSFs become which are always smaller than the original standard deviation.
The relation S22 applies also to the micro-images Notably, the proportionality factor in the above relation is independent of x s . Thus, the above condition is satisfied only if equation S24 can be factored into a term depending only on x s and another term depending only on x d . This can occur only in the following cases: 1. If o is constant over the whole space, namely it is a flat object. This case is trivial and we will not consider it.
. Namely, the shape of the complete PSF h has to be independent of x d . This condition can be satisfied in ideal confocal microscopy.
4. If the illumination spot is so small that only the position x s of the object o contributes to image formation. This case is satisfied only using a pointlike excitation h exc (x d ) = δ(x d ), as in ideal STED microscopy.
We are interested only in the third and fourth case, analyzed in detail in the following sections.

The confocal case
In this section, we assume that the following relation is exact where the complete PSF h[x s − µ(x d )] in general is not normalized. The normalization factor is calculated as Thus, we can rewrite the PSF as where g(x s ) is a normalized function. Using assumption S29, we replicate the calculations that led to equation S23 after performing APR Notably, the scaling factor s(x d ) is exactly the convolution of the excitation PSF with the detection PSF. Thus, we can write the equation of the micro-image after APR as where γ = R 2 o(x s ) dx s is a scale factor that depends only on the sample distribution.
Therefore, we have found that under the assumption S29, the reassigned micro-image is proportional to the fingerprint. Importantly, this mathematical relation is satisfied by a Gaussian PSF, which is well known to be a good approximation for both detection and excitation PSFs if x d 1 Airy Unit. However, our result is more general, being valid for each function that satisfies equation S29. Moreover, APR is an adaptive algorithm which can compensate for mild deviations from ideality, since the shift-vectors are calculated as the translations that maximize the similarity between the scanned images.

The STED case
Assuming an ideal STED system, we can impose a point-like excitation where ε is a (theoretically vanishing) scale factor. In this case, the shift vectors are null Thus, APR does not modify the micro-images. Using assumption S35, equation S23 becomes is a proportionality factor with the units of a photon flux.
We found the interesting result that the micro-image obtained while performing STED microscopy is the detection PSF.
Applying condition S35 to equation S34, we obtain the fingerprint of a STED image which corresponds to the detection PSF, except for a scale factor.

Focus-ISM
We have found that, after APR and under reasonable assumptions, the microimages of an ISM dataset are proportional to the fingerprint of the image. In this section, we generalize the equations for the micro-images and for the fingerprint to the case of a three-dimensional object o(x s , z). The micro-image of a 3D object is where the convolution is calculated along the dimensions (x d , y d ). Assuming that the sample is localized on a single plane, we have o( . Using this assumption, we obtain the micro-image of an object o placed at distance z 0 from the image plane z After pixel reassignment, the micro-image becomes The fingerprint of an image generated by a 3D distribution of emitters o(x s , z) is where γ(z) = R 2 o(x s , z) dx s is a weight function that depends only on the sample distribution, and the first convolution is calculated over the z-axis, while the second convolution is calculated over the x s axes. Using the assumption that the sample is localised on a single plane, we have γ(z) = γ 0 δ(z − z 0 ) and we obtain which is proportional to equation S41. We can use this result to write the equation of a micro-image generated by multiple emitters at different axial planes Notably the fingerprint is always centered with respect to the detector array, but its lateral width depends on the axial position of the emitter. Approximating again each PSF with a Gaussian function, the fingerprint is a Gaussian function as well. Considering only two contributes (in-focus and out-of-focus), we finally obtain where the weights of the Gaussian functions respect the conservation of the energy constraint Notably, the standard deviation of the in-focus component can easily be measured experimentally. The fingerprint is effectively a micro-image averaged over many scan points. Thus, the SNR of the fingerprint is much higher than that of a micro-image. We can exploit this information to calibrate σ sig . This process can be performed with a separate experiment, or alternatively, it is possible to calculate the fingerprint on a region of i(x s |x d ) containing only in-focus emitters. If fitted to a Gaussian function, the fingerprint of the sub-image provides an experimental measure of σ sig . Once each micro-image has been separated in an in-focus and out-of-focus component, the corresponding images are constructed as follows Supplementary References   The signal is measured as the peak value of the images of the beads. The graph reports the average values and the corresponding standard error. Note that at high depletion power the SNR of the raw dataset decreases, potentially hindering the successful ISM reconstruction. Thus, the peak signal level of the ISM images could appear lower than that of the original images.   Figure S6: Distribution of light on the SPAD array. From the PSF simulations presented in Figure 4 we calculated a z-stack of fingerprints. In a, we show the ratio of the intensity averaged over the border elements to the intensity of the central element of the fingerprint. In b, we show the standard deviation resulting from a Gaussian fit to the fingerprints. In c, we show the ratio of the standard deviation of confocal fingerprint to the stadard deviation of the STED fingerprint. In the focal plane, the ratio is 1.32 to be compared with the theoretical expected value of 1.38. The axial units are normalized with respect to the depth-of-field (DoF = λn/NA 2 ). Figure S7: STED and confocal micro-images. Simulated micro-images generated by a single-emitter at different lateral and axial positions. In a we imposed ideal STED conditions (ς = 300). The micro-image corresponds to the detection PSF and contains a non-negligible amount of signal only if the emitter is centered in the xy plane. The lateral confinement occurs also out-of-focus. In b, we simulated the confocal case. The micro-images shift as the emitter moves laterally. In c, we show the micro-images from b after adaptive pixel reassignment. The new micro-images are now all centred and independent from the sample structure, except for a scale factor. In detail, they are proportional to the fingerprint of the corresponding axial plane. For visualization purposes, the simulatons are run with a detector size of one Airy unit. Figure S8: PSF in case of annular excitation. In a, we show the confocal PSFs obtained using an annular illumination. We calculated the PSFs for each detector element and for both the lateral and the axial plane. In b, we show the fingerprints calculated by summing all the scan points of the sections of the PSF at different axial positions. Notably, the intensity profile is flatter than in the case of Gaussian beam excitation. Figure S9: Working principle of the Focus-STED algorithm. In a, we show the core of the analysis: each micro-image is fitted to two Gaussian functions. The micro-image is normalized before the fitting and the weights respect the condition α + β = 1. The width of the in-focus Gaussian function is fixed either by using the theoretical value or by an additional measurement. In b, we show how to experimentally measure σ sig . We select a sub-image containing only in-focus emitters and calculate the corresponding fingerprint. The result is fitted to a single Gaussian function, whose standard deviation is used as σ sig for the pixel-by-pixel fitting of the model presented in a. Figure S10: Result of the Focus-ISM algorithm. The first five reconstructions are obtained by fixing σ bkg . For σ bkg = σ sig the in-focus and out-of-focus photons follow exactly the same distribution and cannot be distinguished by the algorithm. Thus, the light is simply split into two identical images. For increasing σ bkg , the algorithm distinguishes the two contributions, but if σ bkg is too small some misclassification may happen and too much signal is removed. For divergent σ bkg , the shape of the background approaches a constant offset. Notably, the result is still superior of that of the simple subtraction because of the physical constraint imposed by the fitting process. The last reconstruction is obtained letting σ bkg be a fitting parameter, with the only physical constraint of being much larger than σ sig . With this last approach, σ bkg is found adaptively pixel-by-pixel.  Bottom left, image of the background removed. The xz and yz slices are taken from the middle of the volume and are linearly interpolated along the zaxis for visualization purposes. All images were acquired using an STED power of 125 mW measured before the objective lens. Figure S13: Focus-ISM applied to a z-stack of confocal images. Top left, confocal image acquired with a closed pinhole. The SBR is high, but the SNR is low. Top right, STED image acquired with an open pinhole. The SNR is high, but the SBR is low. Bottom right, result of the focus-ISM algorithm applied to the same dataset of the above images. Both the SNR and SBR are high, and the similarity with the closed-pinhole image proves the absence of artifacts. Bottom left, image of the background removed. The xz and yz slices are taken from the middle of the volume and are linearly interpolated along the z-axis for visualization purposes. Figure S14: Radial spectra of experimental images. The spectra have been obtained by calculating the Fourier transform of the XY images shown in Figures 6, S12, and S13, performing the angular average, and normalizing the result with respect to the mean of the image. The legend contains the value of the integral of the spectra over the spatial frequencies k. A higher value of the integral can result from an increased spatial resolution -namely if the cut-off frequency is higher -or if the contrast at high frequencies is enhanced with respect to the low frequencies. a, the superior lateral resolution of ISM, compared to LSM (1.40 AU), determines a larger a cut-off frequency and a larger integral value. Although having poor SNR, LSM (0.28 AU) has superior lateral resolution and greater optical sectioning, thus showing higher contrast of large spatial frequencies and a bigger integral value. Focus-ISM has good SNR, the same lateral resolution of ISM, and the best optical sectioning capability. These features result in the biggest integral value. b, STED (1.40 AU) and STED-ISM have a comparable radial spectrum, due to the high power of the depletion beam. STED (0.28 AU) has better high-frequency contrast, but the best contrast is still achieved by Focus-STED-ISM. Figure S15: Collection efficiency at different pinhole sizes. In a, we show the results of the simulated collection efficiency from an in-focus sample at different pinhole size. In the STED case, the distribution of the light on the detector plane is given by the detection PSF. In the confocal case the light distribution is given by the convolution of the excitation PSF with the detection PSF, thus broader of a factor √ 2 at the focal plane. Thus, the correct estimation of the correct pinhole size should consider the full fingerprint. The detection PSF is a an underestimation of intensity distribution on the detector, except than in the ideal STED case. In b, we show the full width at half maximum (FWHM) of the confocal PSF as a function of the diameter of the pinhole. The resolution and the collected signal are inversely proportional.   Fig. 4. a, 2D MTFs of LSM (pinhole size: 1.40 and 0.28 AU), ISM and focus-ISM, plotted over the axial (k z ) and lateral (k x ) spatial frequencies. The MTF of LSM with an open pinhole (1.40 AU) and the MTF of ISM differ only in the extension along the k x -axis, due to the improved lateral resolution of ISM, but they both show the missing cone along the k z -axis, namely the lack of optical sectioning. The MTF of LSM with a closed pinhole (0.28 AU) and of focus-ISM present the same extension along the k x -axis, and a similar extension along the k z -axis, without the presence of the missing cone. b, 1D Normalized MTFs, obtained by integrating the corresponding 2D MTFs of a along the axial or lateral dimension. The axial profile clearly shows the superior optical sectioning of focus-ISM compared to LSM (0.28 AU), while keeping a signal content comparable to that of LSM (1.40 AU). The lateral profile shows the enhanced lateral resolution of ISM and focus-ISM with respect to LSM (1.40 AU).  Fig. 4. a, 2D MTFs of STED (pinhole size: 1.40 and 0.28 AU), STED-ISM and focus-STED-ISM at saturation factor ς = 300. b, 1D Normalized MTFs, obtained by integrating the corresponding 2D MTFs of a along the axial or lateral dimension. In this case all the lateral MTFs have the same support on the k x -axis because of the high-value of ς. The axial MTFs show that Focus-STED-ISM still improves the optical sectioning, as shown by the higher contrast of the axial MTF compared to the others.