Achromatic light patterning and improved image reconstruction for parallelized RESOLFT nanoscopy

Fluorescence microscopy is rapidly turning into nanoscopy. Among the various nanoscopy methods, the STED/RESOLFT super-resolution family has recently been expanded to image even large fields of view within a few seconds. This advance relies on using light patterns featuring substantial arrays of intensity minima for discerning features by switching their fluorophores between ‘on’ and ‘off’ states of fluorescence. Here we show that splitting the light with a grating and recombining it in the focal plane of the objective lens renders arrays of minima with wavelength-independent periodicity. This colour-independent creation of periodic patterns facilitates coaligned on- and off-switching and readout with combinations chosen from a range of wavelengths. Applying up to three such periodic patterns on the switchable fluorescent proteins Dreiklang and rsCherryRev1.4, we demonstrate highly parallelized, multicolour RESOLFT nanoscopy in living cells for ~100 × 100 μm2 fields of view. Individual keratin filaments were rendered at a FWHM of ~60–80 nm, with effective resolution for the filaments of ~80–100 nm. We discuss the impact of novel image reconstruction algorithms featuring background elimination by spatial bandpass filtering, as well as strategies that incorporate complete image formation models.


D 1
Dichroic 390/485/593 nm bandpass (Chroma ZT390/485/593tpc-UF5) D 2 Dichroic 580 nm edge BrightLine (Semrock FF580-FDi01-25x36) D 3 Dichroic 387 nm long pass (Chroma T387lp-UF1)  (a) An acquired camera image in a small region of 120×90 pixels (about 35×26 illumination periods). The analog-to-digital converter (ADC) bias was subtracted and the ADC conversion factor was applied to scale the image in photo-electrons. (b) Bandpass filtered image obtained by convolving the acquired image a with the Laplacian-of-Gaussian filter kernel shown to scale in-between the images. (c) Determination of the pinholes and sensitivities. Left: Simulated image of equally bright emitters at every position of a null. Right: Bandpass filtered simulated image with applied zero threshold. Positive values are identifying the pinhole regions (4-6 pixels within a 3×3 neighbourhood at the nulls). (d) Extracted signals at the nulls obtained by the integrated signal in each pinhole normalized with its sensitivity. (e) The pinhole PSFs are extracted from simulated and analysed images. Equally bright emitters are placed at every fifth null to build all PSFs with only 5×5 images. (f) The PSF matrix is then applied in a least squares conjugate gradient fit to find the closest matching pinhole signals (shown) that yield the observed image d. . A sequence of two images each, of dense keratin19-Dreiklang structures, were taken, with activation, off-switching and read-out illuminations. For the first image, the pattern was displaced by half a period (180 nm) between the off-switching and the read-out step, as would be suitable if these patterns overlapped completely. For the second image, an additional variable shift of 10-360 nm in steps of 10 nm was introduced. Fluorescence in the recorded images was therefore modulated due to variable spatial coincidence of switch-off and read-out, with it being minimal when the two patterns (off-switching and read-out) were completely overlapping, and being maximal when completely out-of-phase. From this ratio we determined the relative offset between the patterns of two wavelengths and took it into account for later measurements. Black dots represent the ratios of the mean fluorescence in the images with the offset to the mean fluorescence in the images without the offset. Red lines represent the fit to sinusoidal functions with 360 nm period. The maximal ratio was found at the optimal offsets of 320 nm (a) and 280 nm (b).

Optical setup
We adapted the wide-field RESOLFT instrument by Abberior Instruments, whose original optical setup is largely based on setup described in [1]. The modified RESOLFT microscope is shown in Supplementary Figure 1. An Olympus IX-83 microscope body was used and customized at its side port. Three continuous-wave laser beams of 405 nm, 488 nm and 592 nm wavelengths were combined by clean-up filters and dichroic mirrors and coupled into a polarization-maintaining single-mode fiber. An achromatic lens L 1 collimated the beams at the fiber output to a beam diameter of about 12 mm. An achromatic half-wave retarder plate (λ/2) and polarizing beam splitters (PBS) sent the beams at equal powers through two binary phase gratings (G 1 and G 2 ). The phase gratings diffracted the spolarized beams primarily into the first diffraction orders. These were focused by the achromatic lens L 2 and selected by an order selection mask (OS) to cancel stray light and the order zero. The achromatic lenses L 3 and L 4 relayed the illumination beams on a tip-tilt piezo scan mirror. The custom-made dichroic mirror D 1 transmitted the wavelength ranges of 350-430 nm, 480-490 nm and 588-598 nm. The achromatic lens L 5 and the microscope tube lens relayed the beams into the entrance pupil of the microscope objective. Thereby, 50× demagnified images of the phase gratings were produced at the objective's focal plane. Hence, the sample was illuminated with square grids of 360 nm periods. When imaging Dreiklang, we switched the protein's fluorescence on with a homogeneous illumination of 365 nm wavelength provided by a LED source and injected via the dichroic mirror D 3 . The fluorescence was collected by the microscope objective and sent back via the scan mirror to the main dichroic mirror D 1 . The lenses L 5 and L 4 and the achromatic lenses L 6 and L 7 /L 8 relayed the intermediate image onto the sensors of two scientific CMOS cameras. The wavelength ranges of 440-470 nm, 500-575 nm and 610-750 nm were reflected by D 1 and further split by the dichroic mirror D 2 into two colour channels. The dichroic mirrors and the bandpass filters F 1 and F 2 selected the fluorescence at 500-575 nm and 610-680 nm wavelengths.

RESOLFT image reconstruction
The raw wide-field images were analysed in three major steps: 1. Determination of the signal origins corresponding to the positions of the intensity minima ("nulls") of the switching off illumination pattern. 2. Estimation of the in-focus fluorescence emission from the regions of the nulls. 3. Reattribution of the fluorescence emission to confined regions at the nulls.
We improved the image reconstruction over the previously published method [1] to (i) more accurately identify the positions of the nulls by accommodating distortions of the illumination pattern and/or its images, and to (ii) better reject the fluorescence emission from out-of-focus features. The analysis can exclude sensor pixels with high dark current or excessive read noise as identified in dark pictures. The new analysis algorithms are detailed below.

Signal origins
Whenever possible, the set of positions , of the nulls were extracted from the raw images of a measurement (Supplementary Figure 2). For samples with very sparse features, an extra calibration measurement of a thin layer of fluorescent proteins was performed shortly prior to the image acquisition of the sample.
We assume that the nulls are located on an approximately regular square grid aligned approximately with the camera sensor along its horizontal and vertical directions. The nulls are projected at the same sensor positions for the image frames captured during a measurement. Therefore, the raw images of an entire acquisition were summed to obtain , , an image of the average fluorescence emission in all unit cells across the field of view (Supplementary Figure 2a).
The grid periods Λ and Λ were estimated by identifying the spatial frequency of the major peak in the horizontal and vertical spatial power spectra of the sum image (Supplementary Figure 2b). These spectra were calculated with the Fourier transform with a frequency sampling fine enough to achieve a precision of about three digits in the period estimates (we used 128× oversampling). To account for potential misalignments and distortions of the illumination grid, the power spectra were calculated on moving averages spanning 30 pixels along the opposite direction (Supplementary Figure 2a). The periods were then estimated by a weighted average on the peak spatial frequencies, where the weights were proportional to the peak spectral power density (Supplementary Figure 2b).

Signal and background estimation by local pinholes
Signals were locally integrated using Gaussian-weighted masks to include neighbouring pixels at the nulls. For each scan step, the local background at each null was determined by 80% of the mean of the signals obtained at the four nearest-neighbour positions diagonally in between the nulls, and then subtracted. For further details of the computational image reconstruction see [1].

Background elimination by spatial bandpass filtering
Supplementary Figure 3 illustrates the image analysis. Each raw image (Supplementary Figure 3a) was smoothed with a Gaussian filter of FWHM Λ Λ / 4√2 for reducing noise and the Laplacian was taken from the smoothed images to eliminate the background (Supplementary Figure 3b). These bandpass-filtered raw images were then integrated in small regions (pinholes) at the nulls' positions to estimate the fluorescence emitted there (Supplementary Figure 3d). The pinholes were determined by the positive signals in 3×3 camera pixels neighbourhoods centred at the nulls (Supplementary Figure 3c). In order to equalize their sensitivities, the pinholes' signals were determined by processing a simulated image of equally bright emitters at the nulls.
This Laplacian-of-Gaussian filtering eliminated background efficiently but introduced some artefacts due to partially overlapping images of adjacent nulls (observe for instance negative values in Supplementary Figure 3b,d). Therefore, images of fluorescence emitters in the focal section at single nulls were simulated (described in the next section) and processed to estimate the point spread functions including the image analysis (Supplementary Figure  3e). The artefacts were then reduced by redistributing the pinhole signals (Supplementary Figure 3c) such that the squared difference to the observed bandpass filtered image (Supplementary Figure 3b) was minimized (Supplementary Figure 3f).

Model-based maximum likelihood estimation of signals
Given the parameters of the microscope objective and the illumination optics, the illumination patterns in a lateral plane of the sample were estimated with the fast focus field calculation tool [3] (Supplementary Figure 4).
The on-and off-switching and the emission of the fluorescent protein is illustrated in Supplementary Figure 4a and was estimated as described in the RESOLFT section in reference [4]. We assumed to switch on up to 80% of the proteins ( 1.6 . 1 We assumed an off-switching illumination dose to reach a five-fold improvement in lateral resolution ( 25), which corresponded to 60-70 nm FWHM of the regions with emissive fluorescent proteins. We further assumed to switch off up to 80% during read-out of the fluorescence ( 1.6) for rsCherryRev1.4 and up to 25% ( 0.3) for Dreiklang. 2 The detection point spread function (PSF) was estimated with the calculation described by Leutenegger and Lasser [5]. The detection PSF was first calculated in lateral planes with a resolution of 2 nm over an axial range of 2 μm in steps of 10 nm. Neighbouring lateral planes were then grouped into sections such that the areas of the PSFs in adjacent sections quadrupled (Supplementary Figure 4b). The planes of the focal section were grouped such that the focal PSF area increased by 42% with respect to the smallest PSF area in the focal plane. Thus, a set of significantly different PSF sections PSF , was obtained, centered at axial positions and spanning an axial range Δ . In particular, we used ∈ 110,305,680,1390 nm and Δ ∈ 250,120,650,780 nm for modelling the focal section at the glass cover slip and three neighbouring sections in the sample. The fluorescence signal from a unit cellcalculated as outlined in the previous paragraph -was then convoluted with the section PSFs to estimate the images of features in these sections.
The image formation model consisted of a sparse matrix , whose element was the contribution to pixel by the PSF of the coefficient. The elements of the coefficient column vector were the estimated fluorescence signals at the positions in the sections . The measured image in column vector was thus obtained by the matrix product plus shot noise and read-out noise . In order to build the model column by column, the section PSFs were placed at the lateral positions and binned into the camera pixels. Thereby, subpixel shifts of the nulls' positions with respect to the camera pixels were taken into account.
A conjugate gradient solver was used to find the coefficients maximizing the Poisson likelihood (photo-electron statistics) that the estimated image corresponds to the measured raw image . This fitting procedure eliminated the blurred signals of defocused features. However, only the coefficients of the focal section could be estimated with high spatial resolution as the others showed significant artefacts due to the overlap of the signals of adjacent nulls. The iterative fitting took about an hour per 100 raw images of a RESOLFT image acquisition.

Signal reattribution
The retrieved signals were combined into the RESOLFT image by placing each signal at its new pixel position in the final image. Photo-bleaching was accounted for by weighing the signals with the inverse of the average signals in each raw image, whereby a second-order polynomial decay of the average signals was assumed to avoid artefacts due to the sample structure.

Wide-field image reconstruction
In order to compare the RESOLFT images with wide-field images, we combined the raw images as follows. We upsampled the raw images to match the spatial sampling of the RESOLFT image and shifted each raw image to its position within the unit cell of the illumination patterns. These steps were performed by zero-padding and phaseshifting the spatial spectra of the raw images. The resampled raw images were then summed with weights inversely proportional to the average signals in each image to account for photo-bleaching (see 3. Signal reattribution). This wide-field image was further processed by a Lucy-Richardson deconvolution, where we used the theoretical diffraction-limited detection efficiency as PSF.