Image filtering in structured illumination microscopy using the Lukosz bound

Christiaan H. Righolt; Johan A. Slotman; Ian T. Young; Sabine Mai; Lucas J. van Vliet; Sjoerd Stallinga

doi:10.1364/OE.21.024431

1. Introduction

The diffraction limit to resolution is broken in optical microscopy by numerous techniques developed in recent years [1, 2]. Localization microscopy (PALM, STORM) circumvents the diffraction limit by localizing switchable fluorophores with nanometer precision, STED improves resolution by reducing the fluorescence emission to a sub-diffraction sized region by de-exciting fluorophores near the boundary of the excitation spot with a ring shaped depletion spot. The technique of Structured Illumination Microscopy (SIM) extends the diffraction limit by a factor of two by illuminating the sample with a series of periodic patterns of closely spaced lines of different orientations and phases. The final super-resolution image is reconstructed from this set of recorded images.

The idea to use heterodyne detection to enhance the effective spatial resolution in microscopic imaging was introduced in the 1960s by Lukosz [3]. The application of these techniques to achieve optical sectioning rather than resolution improvement were further developed in the mid-1990s [4, 5]. Implementations to increase the effective lateral resolution were developed around the turn of the century [6–9]. Two-dimensional SIM was further extended by exploiting the non-linear fluorescence response [10, 11] and by combination with solid immersion imaging [12]. The principles of SIM were implemented in 3D several years ago [13]. More recently the temporal resolution of the technology has improved in both two- [14] and three- [15] dimensional SIM, based on the use of Spatial Light Modulators (SLMs) [16]. A reduction of out-of-focus light can be achieved using a combination with line scanning [17].

Image reconstruction in SIM boils down to first disentangling the different spatial frequency bands, and subsequently stitching these bands together to extend the bandwidth of the reconstructed image. This step is often integrated with a filtering approach for suppressing noise and for tailoring the effective response as a function of spatial frequency. Commonly used filters in SIM image reconstruction are borrowed from the field of image restoration, i.e. methods to extend resolution beyond the diffraction limit by incorporating prior knowledge such as nonnegativity [18], noise models [19], background level [20], and image smoothness [21]. The default filtering technique in SIM is (generalized) Wiener filtering [22–24], with a signal-to-noise ratio (SNR) that is assumed constant across the entire spatial frequency range. The free parameter in this type of filter is the so-called regularization parameter, which measures the ‘strength’ of the prior knowledge. It is not a priori clear what value this parameter should have. Moreover, if the parameter takes a small value, several artefacts appear, in particular edge ringing, halos surrounding bright objects, and negative pixel values. The origin of the artefacts is an effective transfer function which is flat up to a point close to the extended diffraction limit. A suitable counter measure is therefore the application of a so-called apodization filter, which suppresses the response at high spatial frequencies compared to the response at low spatial frequencies [8, 13]. Apodization filters are applied ad-hoc and often introduce further free parameters which have optimal settings that are a priori unclear. An additional image artefact common to SIM is the enhancement of low-frequency noise structures in the reconstructed image. The effect of the different free parameters on this artefact have not been explored so far.

The goal of the current paper is to investigate the role of filtering and filter parameters on the image quality of SIM. The issues of edge ringing and negative Point Spread Function (PSF) values were treated systematically by Lukosz half a century ago [25, 26]. He derived a set of conditions on the Optical Transfer Function (OTF) of the system that must be satisfied in order to guarantee a non-negative PSF. This Lukosz bound to the OTF thus provides a necessary (but not a sufficient) condition for having a non-negative PSF. Here, we propose to use the Lukosz bound for apodization in SIM image reconstruction and filtering. We will present and evaluate different ways to incorporate the Lukosz bound. The advantage of our proposal is that it provides a systematic approach as opposed to ad-hoc procedures and that it explicitly rules out free apodization parameters.

The content of the paper is as follows. The key results of Lukosz are summarized in section 2, image formation and filtering for SIM is presented in section 3. The application of the Lukosz bound in image reconstruction is described in section 4. Simulations are presented in section 5, experiments in section 6. The paper concludes with a discussion of the results and an outlook to further work in section 7.

2. The Lukosz bound

Lukosz considered linear shift-invariant band-limited imaging systems, characterized by a PSF H (u), or equivalently its Fourier transform, the OTF Ĥ (v). Here, and thoughout the paper we use the designation u for the spatial coordinates, and v for the spatial frequencies. The PSF and OTF are a Fourier transform pair, for which we use the convention:

\hat{H} (v) = \int H (u) exp (- 2 π i u \cdot v) d u .

We will first consider 1D band-limited imaging systems, with a cut-off spatial frequency q_c. Lukosz showed in his papers [25,26] that the Modulation Transfer Function (MTF), the absolute value of the OTF, must be smaller than the Lukosz bound in order to have a non-negative PSF, i.e. a PSF that is positive-definite everywhere:

| \hat{H} (v) | \leq Λ_{1}^{stairs} (v) = cos (\frac{π}{M + 1}) for | v | \geq \frac{q_{c}}{M} with M = 1, 2, 3, \dots

This function describes a staircase with infinite, progressively smaller, steps toward the origin, as shown in Fig. 1. It is a necessary condition that the MTF should be below the ‘Lukosz stairs’, but not a sufficient condition. An additional condition is that the average of the OTF over its support must be less than 1/2. A suitable curve in practice is the continuous curve connecting the lower corners of the stairs:

{\hat{Λ}}_{1} (v / q_{c}) = cos (\frac{π | v | / q_{c}}{| v | / q_{c} + 1}) for | v | \leq q_{c} and {\hat{Λ}}_{1} (v) = 0 elsewhere,

This practical Lukosz bound OTF is shown in Fig. 1, as well as the incoherent OTF and the corresponding Lukosz bound PSF and incoherent PSF. It turns out that indeed a positive definite PSF is obtained, at least for spatial frequencies up to the second diffraction minimum. The achieved sharpness is only marginally less than the sharpness of the incoherent PSF.

Fig. 1 (a) Plot of the Lukosz bound OTF and incoherent OTF in 1D. (b) Plot of the Lukosz bound PSF and incoherent PSF in 1D. (c) Plot of the Lukosz bound OTF and incoherent OTF in 2D. (d) Plot of the Lukosz bound PSF and incoherent PSF in 2D.

Download Full Size | PDF

Next, consider a 2D imaging system. Now, the cut-off spatial frequency q_c (ϕ) can depend on the azimuthal angle ϕ. According to Lukosz the MTF for a spatial frequency vector (v cosψ, v sinψ) should be below the product of the two 1D Lukosz bounds along the orthogonal axes (cosϕ, sinϕ) and (−sinϕ, cosϕ) for any value of ϕ. The 2D bound for azimuthal angle ϕ is thus:

{\hat{Λ}}_{2}^{azi} (v, ϕ) = {\hat{Λ}}_{1} (v_{x^{'}} / q_{c} (ϕ)) {\hat{Λ}}_{1} (v_{y^{'}} / q_{c} (ϕ + π / 2)),

where the spatial frequencies in the rotated coordinate frame are v′ = (v cos(ψ −ϕ), v sin (ψ − ϕ)). The MTF should be below this bound for any value of ϕ giving the overall 2D Lukosz bound as:

{\hat{Λ}}_{2} (v) = min {{\hat{Λ}}_{2}^{azi} (v, ϕ) | ϕ \in [0, 2 π)} .

The cut-off spatial frequency is constant as a function of azimuthal angle in case the OTF is rotationally symmetric. In that case the recipe for the 2D Lukosz bound reduces to:

{\hat{Λ}}_{2} (v) = min {cos (\frac{π | v |}{| v | + q_{c}}), {cos}^{2} (\frac{π | v |}{| v | + \sqrt{2} q_{c}})} for | v | \leq q_{c}

and 0 elsewhere, i.e. the minimum is taken between the product of the two 1D-Lukosz bounds for coordinate frames rotated over angles 0 and π/4. Fig. 1 shows the 2D OTF and PSF according to the Lukosz bound of the rotationally symmetric imaging system, as well as the conventional incoherent OTF and PSF. The Lukosz bound PSF is positive-definite up to the second dark diffraction ring (where it drops slightly below zero to a value around −0.006) and is a bit sharper than the incoherent PSF (full-width at half-maximum (FWHM) is around 5% smaller). The possibly negative pixel values arising from the small negative PSF values around the second dark diffraction ring are not expected to give rise to significant edge ringing or halo artefacts. Also, the effects are mitigated by a background signal which is always non-zero in practice and by blurring through the finite pixel size.

In case the OTF of an imaging system is equal to the Lukosz bound, it is intuitively clear that this provides the most compact non-negative PSF that can be achieved. In a sense the Lukosz bound, therefore, describes the ‘ideal’ imaging system. It is further noted that the Lukosz bound does not improve much over the standard incoherent response.

3. Image formation in SIM

Image formation in SIM consists of two distinct steps, both of which will be discussed in this section; 1) the actual physical image acquisition through the microscope, and 2) the postprocessing that yields the final reconstructed image. Here, the spatial coordinates in object space are normalized with λ/NA_ob with NA_ob the objective numerical aperture, while the coordinates in image space are normalized with λ/NA_tube with NA_tube the tube lens numerical aperture. The imaging system has unit magnification in these normalized coordinates. The spatial frequencies are normalized accordingly with factors NA/λ. The stop of the imaging system thus corresponds to the unit circle, and the normalized cut-off spatial frequency of a conventional widefield fluorescence microscope is then equal to two. The following functions will be used throughout: T (u), T̂ (v) for the object function/spectrum; W (u), Ŵ (v) for the illumination function/spectrum; and H (u), Ĥ (v) for the PSF/OTF of the widefield fluorescence microscope. The symbol ⊗ will denote convolution.

3.1. Image acquisition and phase shifting

The premise of SIM is that a periodic illumination function is used that is rotated and shifted with respect to the sample. In particular, a periodic pattern of stripes is used:

W (u) = W_{av} [1 + 2 w cos (4 π q u_{x})],

with W_av the average illumination intensity, w a parameter characterizing the modulation depth (optimum value w = 1/2), and where p = 1/(2q) is the spacing of the lines. Such an illumination pattern can be made from the interference pattern of two plane waves traveling at opposite angles with respect to the optical axis. The illumination pattern in the back focal plane of the objective lens (assuming epi-illumination is used) then consists of two dots at positions R(q, 0) and R(−q, 0), with R the pupil radius, from which it follows that q ≤ 1. The Fourier transform of the illumination function consists of a set of discrete delta-peaks:

\hat{W} (v) = \sum_{m} w_{m} δ (v - q_{m}),

with q_m = m(q, 0) and with non-zero Fourier coefficients (w₋₂, w₀, w₂) = W_av (w, 1, w). A set of images is recorded for rotations R_l, (l = 1, 2, . . . N_r) and translations u_n = n(p/N_t, 0) (n = 1, 2, . . . N_t). Typically N_r = 3 rotations and N_t = 3 translations are used (for 2D-SIM). The recorded image for rotation l and translation n is given by:

I_{rec, \ln} (u) = H (u) \otimes [W (R_{l} \cdot u + u_{n}) T (u)] .

This can be rewritten in the spatial frequency domain as:

{\hat{I}}_{rec, \ln} (v) = \hat{H} (v) \sum_{m} Φ_{n m} w_{m} \hat{T} (v - R_{l} \cdot q_{m}),

with phase shifting matrix:

Φ_{n m} = exp (2 π i q_{m} \cdot u_{n}) .

The Fourier transform of each recorded image is a linear sum of bands in spatial frequency space. These bands can be distilled by an inversion procedure:

{\hat{I}}_{bands, lm} (v) = \sum_{n} Φ_{n m}^{*} {\hat{I}}_{rec, \ln} (v) = w_{m} \hat{H} (v) \hat{T} (v - R_{l} \cdot q_{m}) .

This inversion works provided the number of translations N_t ≥ K, with K the number of independent Fourier components of the illumination function. In case the translations are not equidistant, the Moore-Penrose pseudo-inverse can be used to obtain the band images.

3.2. Weighted sum of shifted bands reconstruction

The most simple and straightforward method to obtain a final super-resolution image is to take a weighted linear combination of the bands shifted in spatial frequency space [6, 7, 27, 28]:

{\hat{I}}_{lin} (v) = \sum_{lm} s_{m} {\hat{I}}_{bands, lm} (v + R_{l} \cdot q_{m}) = \sum_{lm} w_{m} s_{m} \hat{H} (v + R_{l} \cdot q_{m}) \hat{T} (v) .

Here (s₋₂, s₀, s₂) = S_av (s, 1, s) are a set of weight coefficients, which can in principle be chosen at will. The parameter S_av can be used to scale the overall effective OTF:

{\hat{H}}_{lin} (v) = \sum_{lm} w_{m} s_{m} \hat{H} (v + R_{l} \cdot q_{m}),

such that Ĥ_lin(0) = 1. The parameter s can be used to tune the relative height of the side bands compared to the central band. It is common to pick this side band height parameter a priori - often as s = 1 or s = 1/w - and not estimate it as part of a filtering process [8, 13, 24]. Interestingly, the band stitching process can be expressed as a convolution operation in spatial frequency space. It follows that this is equivalent to point multiplication of the recorded images with a suitable mask function in real space, and then adding all contributions:

I_{lin} (u) = \sum_{\ln} S (R_{l}^{- 1} \cdot u + u_{n}) I_{rec, \ln} (u),

with mask function:

S (u) = S_{av} [1 + 2 s cos (4 π q u_{x})] .

So, the mask function for each recorded image has the same period, orientation and phase as the illumination pattern for the recorded image, only the modulation depth is characterized by the parameter s instead of w.

3.3. Filtering approaches

A next step in sophistication is to apply a noise suppression filter with kernel F̂ (v) to the image obtained by the weighted sum of shifted bands reconstruction method, producing a super-resolution image Î_apo (v) = F̂ (v)Î_lin (v). It urns out that the level of noise suppression for this post filtering approach can be improved by integrating the filtering operation into the process of stitching together the bands in spatial frequency space. This is done by making the filter kernel dependent on the spatial frequency band. This so-called generalized filtering approach gives rise to:

{\hat{I}}_{gen} (v) = \sum_{lm} s_{m} {\hat{F}}_{lm} (v) {\hat{I}}_{band, lm} (v + R_{l} \cdot q_{m}) .

The filter kernels can be derived from minimizing a Tikohonov-Miller (TM) type of functional [29]:

ℒ_{gen} = \int [\sum_{lm} {| {\hat{B}}_{lm} (v) {\hat{I}}_{gen} (v) - s_{m} {\hat{I}}_{band, lm} (v + R_{l} \cdot q_{m}) |}^{2} + κ {| \hat{A} (v) {\hat{I}}_{gen} (v) |}^{2}] d v .

where the filtered super-resolution image Î_gen (v) is considered as the variable in the minimization procedure. The first term of the TM-functional is the ‘data misfit term’, i.e. the summed squared difference between the data and the forward model, and the second term is the ‘regularization term’, representing the prior knowledge of the reconstructed object. Here the functions B̂_lm (v) describe the ideal or desired response of the imaging system. They are usually identified with the weighted and shifted microscope OTF s_mw_mĤ (v + R_l · q_m), but can in principle be chosen differently. The parameter κ is the regularization parameter and Â (v) is the regularization weight function. We consider:

\hat{A} (v) = {| v |}^{p},

with p an integer, so for p = 0 the total signal energy is regularized, for p = 1 the total signal gradient energy, for p = 2 the total signal Laplacian energy, etc. Minimization of the TM-functional gives:

{\hat{F}}_{lm} (v) = \frac{{\hat{B}}_{lm} {(v)}^{*}}{κ \hat{A} {(v)}^{2} + \sum_{lm} {| {\hat{B}}_{lm} (v) |}^{2}},

and an OTF:

{\hat{H}}_{gen} (v) = \frac{\sum_{lm} s_{m} w_{m} {\hat{B}}_{lm} {(v)}^{*} \hat{H} (v + R_{l} \cdot q_{m})}{κ \hat{A} {(v)}^{2} + \sum_{lm} {| {\hat{B}}_{lm} (v) |}^{2}},

The conventional choice for the functions B̂_lm (v) = s_mw_mĤ (v + R_l · q_m) gives rise to an OTF that is flat up to the cut-off in the limit of a small regularization parameter κ. Such OTF’s do not satisfy the Lukosz bound and thus give rise to edge ringing, halo artefacts and negative pixel values.

We point out that these TM-filters correspond to (generalized) Wiener-filters if the signal-to-noise power spectrum is used as the regularization function [30]. In the literature on SIM these filters with total signal energy regularization (p = 0) are, however, usually referred to as (generalized) Wiener-filters [13]. We will follow this nomenclature, although we believe it is basically incorrect.

3.4. Data misfit weight function

Light originating from out-of-focus layers gives rise to a blurred background and hence a reduced contrast in the captured images. Moreover it may lead to honeycomb artefacts in the reconstruction because the low frequency blur is transfered to the regions close to the side band peaks in spatial frequency space. Suppression of out-of-focus light can be done by subtracting a low-pass filtered version of the acquired images before feeding the acquired images into the image reconstruction procedure. An alternative is to modify the filters with a weight function for suppressing the effect of the data misfit function on the final filter close to the side band peaks in spatial frequency space [31]. It follows that the regularization dominates over the data misfit term for the spatial frequencies for which the out-of-focus light has a significant influence. Here we show that this modification of the filtering step can be incorporated into the framework of TM-reconstruction. The data misfit term in the TM-functional must be changed in such a way that the relative weight of that term compared to the regularization term decreases with decreasing spatial frequency. The TM-functional for generalized filtering then becomes:

ℒ_{gen} = \int [\sum_{lm} \hat{g} (v + R_{l} \cdot q_{m}) {| {\hat{B}}_{lm} (v) {\hat{I}}_{gen} (v) - s_{m} {\hat{I}}_{band, lm} (v + R_{l} \cdot q_{m}) |}^{2} + κ {| \hat{A} (v) {\hat{I}}_{gen} (v) |}^{2}] d v,

with ĝ (v) a function that increases with increasing spatial frequency. For example, a Gaussian data misfit weight gives rise to:

\hat{g} (v) = 1 - α exp (- {| v |}^{2} / 2 σ^{2}),

where the magnitude α and width σ are in principle adjustable parameters. Minimization of the TM-functional leads to:

{\hat{F}}_{lm} (v) = \frac{\hat{g} (v + R_{l} \cdot q_{m}) {\hat{B}}_{lm} {(v)}^{*}}{κ \hat{A} {(v)}^{2} + \sum_{lm} \hat{g} (v + R_{l} \cdot q_{m}) {| {\hat{B}}_{lm} (v) |}^{2}},

and an OTF:

{\hat{H}}_{gen} (v) = \frac{\sum_{lm} s_{m} w_{m} \hat{g} (v + R_{l} \cdot q_{m}) {\hat{B}}_{lm} {(v)}^{*} \hat{H} (v + R_{l} \cdot q_{m})}{κ \hat{A} {(v)}^{2} + \sum_{lm} \hat{g} (v + R_{l} \cdot q_{m}) {| {\hat{B}}_{lm} (v) |}^{2}},

A consequence of relatively large values of the Gaussian magnitude α (close to one) is that dips in the effective OTF arise at the spatial frequencies around the center of the main band and of the side bands. The dips can thus lead to violations of the Lukosz bound, as that is a monotonically decreasing function. Although the resulting OTF does not lead to an imaging system that gives rise to a non-negative output for any conceivable object structure, it may be expected to do so for the particular case where a large out-of-focus contribution is present, i.e. for thick samples. The Gaussian function should, therefore, be tailored to the thickness and density of the sample in question by varying α and/or σ for each case.

4. The Lukosz bound in SIM filtering

The filtering approaches presented so far are clearly inappropriate and some form of apodization is in order. According to Lukosz the undesired artefacts are eliminated if the final OTF is below the Lukosz bound. In this section, we will describe how to incorporate the Lukosz bound in the filtering step.

4.1. 2D SIM Lukosz bound

The bound Λ̂ (v) for (2D) SIM can be derived along the lines presented in sect. 2 provided the spatial frequency cut-off is known. For the case of SIM this corresponds to the combined support of the central band and the side bands. For the typical case of N_r = 3 rotations this results in a ‘flower’-shaped support parametrized by:

q_{c} (ϕ) = 2 q cos ϕ_{r} + 2 \sqrt{1 - q^{2} {sin}^{2} ϕ_{r}},

where q ≤ 1 determines the pitch of the illumination pattern and where:

ϕ_{r} = mod {ϕ + \frac{π}{6}, \frac{π}{3}} - \frac{π}{6}

is the reduced azimuth angle corresponding to each circle segment of the flower-shaped support. Fig. 2 shows the OTF and PSF according to the 2D-SIM Lukosz bound for a value q = 0.9. The support of the OTF is reduced somewhat to a hexagonal shape, the spatial frequency cut-off in the x-direction is 2(1 + q), the spatial frequency cut-off in the y-direction is

\sqrt{3} (1 + q)

. The PSF resulting from the Lukosz bound does not differ much from the incoherent PSF for a spatial frequency cut-off 2(1 + q) with a FWHM of 0.28 (in units λ/NA) in both the x and y-direction. The first dark diffraction ring has a hexagonal imprint with 6 ‘pits’, where the PSF drops slightly below zero (minimum value −1 × 10⁻²). This small violation of positive definiteness is not expected to give rise to highly visible artefacts. However, the hexagonal imprint may be reduced by using an anisotropically stretched version of the Lukosz bound:

\hat{Λ} {(v)}^{'} = \hat{Λ} {(v)}^{(1 + μ cos (6 ϕ_{r})) / (1 - μ)},

where μ is a stretching parameter. The anisotropically stretched Lukosz bound OTF gives rise to a non-negative and isotropic PSF (at least up to the second dark diffraction ring) for a value μ = 0.0656. In the following we will use the non-stretched version of the Lukosz bound for the sake of simplicity.

Fig. 2 (a) Plot of the 2D SIM Lukosz bound OTF. (b) Plot of the 2D SIM Lukosz bound PSF. (c) Cross-sections of the 2D SIM Lukosz bound OTF along the x and y-direction, and the incoherent OTF for cutoff 2 + 2q as reference. (d) Cross-sections of the 2D SIM Lukosz bound PSF along the x and y-direction, and the incoherent PSF for cut-off 2 + 2q as reference.

Download Full Size | PDF

4.2. Lukosz bound filtering

There are two ways in which the Lukosz bound can be used in apodization. The first way is simply to apply it as a post-processing apodization filter Λ̂ (v) after a super-resolution image has been generated with the generalized filtering approach. This results in an OTF:

{\hat{H}}_{gen} (v) = \frac{\hat{Λ} (v) \sum_{lm} s_{m}^{2} w_{m}^{2} {| \hat{H} (v + R_{l} \cdot q_{m}) |}^{2}}{κ \hat{A} {(v)}^{2} + \sum_{lm} s_{m}^{2} w_{m}^{2} {| \hat{H} (v + R_{l} \cdot q_{m}) |}^{2}} .

The second way to use the Lukosz bound is to use modified versions of the functions that characterize the ideal or desired response of the imaging system. If we use B̂_lm (v) = s_mw_mĤ (v + R_l · q_m) /Λ̂ (v) then the OTF becomes:

{\hat{H}}_{gen} (v) = \frac{\hat{Λ} (v) \sum_{lm} s_{m}^{2} w_{m}^{2} {| \hat{H} (v + R_{l} \cdot q_{m}) |}^{2}}{κ \hat{A} {(v)}^{2} \hat{Λ} {(v)}^{2} + \sum_{lm} s_{m}^{2} w_{m}^{2} {| \hat{H} (v + R_{l} \cdot q_{m}) |}^{2}} .

Clearly, this alternative route of incorporating the Lukosz bound gives rise to an OTF that is very similar to the one obtained by applying the Lukosz bound as apodization filter after the filtering operation. The difference is that now the regularization term is effectively damped with a factor Λ̂ (v)² resulting in a an effectively smaller regularization effect at higher spatial frequencies. In the limit κ → 0 the effective OTF approaches the Lukosz bound, just as for the generalized Wiener filtering with Lukosz apodization.

This integrated Lukosz filtering approach gives rise to an effective OTF that satisfies the Lukosz bound for all values of the regularization parameter, with the exception of the signal energy regularization case (p = 0). See the Appendix for more details about ensuring the Lukosz bound is not violated in this case.

5. Simulation results

The reconstruction methods presented in the previous sections of this paper are demonstrated here using simulations. We simulated a virtual resolution chart consisting of 512×512 pixels with pixel size 1/16 (in normalized units of λ/NA). The virtual resolution chart, shown in Fig. 3(a) contains various structures (points, crossing lines, bar patterns, uniform blocks) that can be used to judge image quality. Several groups of bar patterns can be identified, the number indicates the pitch in units of λ/16NA. Pattern “4” is thus at the cut-off for widefield imaging, pattern “2” near the cut-off for SIM. The final images were binned 2×2 giving a 256×256 image with pixel size corresponding to Nyquist sampling at the maximum extended cut-off spatial frequency of SIM (equal to 4 in normalized units of NA/λ). The simulated recorded images were corrupted by shot noise, which is taken to be the dominant source of noise. The reconstructions were performed with an implementation of Eqs. (17, 20) in which either the Lukosz bound was used as apodization filter Eq. (29) or in which the Lukosz bound was incorporated in the Tikhonov-Miller functional Eq. (30). All simulations were performed in MATLAB (Math-Works, Natick, USA) using the DIPimage toolbox [32].

Fig. 3 (a) Digital resolution chart used in the simulations. The numbers with the bar patterns indicate the pitch in pixel units, i.e. in units of λ/16NA. (b) Simulated widefield image, (c) Screenshot of a series of SIM image reconstructions with different side band height parameters using a weighted sum of bands method ( Media 1). The screenshot is the reconstruction with high side band height parameter s = 3 showing noise amplification. (d) Screenshot of a series of SIM image reconstructions with generalized Wiener filtering without apodization for different regularization parameter values ( Media 2). Shown is the reconstruction at low regularization parameter κ = 10⁻⁴, which demonstrates halo artefacts and low-frequency noise structures.

Download Full Size | PDF

Fig. 3(a) shows an image of the simulated resolution chart, Fig. 3(b) of the simulated wide-field image where on average the intensity corresponded to 800 photons per (λ/8NA)² area. Key SIM image reconstruction artefacts are shown in Fig. 3(c) ( Media 1), which shows the effect of the side band height parameter s on the results obtained by the weighted sum of shifted bands reconstruction, and Fig. 3(d) ( Media 2), which shows the effect of the regularization parameter κ on the results obtained by the generalized Wiener filtering in the image reconstruction without any apodization. The simple and straightforward weighted sum of shifted bands reconstruction shows resolution improvement at the expense of noise enhancement for increasing values of the side band height parameter. The optimum value for s seems to be in the range 1.0–1.5. Generalized Wiener filtering without apodization results in grave image artefacts. Significant effects of negative pixel values, edge ringing and an enhanced noise structure overlaying the genuine image can be observed, especially for low values of the regularization parameter.

Fig. 4 shows movies of simulated resolution chart images as a function of the regularization parameter κ for generalized Wiener filtering with Lukosz bound apodization and for integrated Lukosz filtering with p = 0 and p = 1 regularization for a side band height parameter s = 1. Both Lukosz bound apodization and integrated Lukosz filtering largely solve the aforementioned image artefacts. Notice that residual noise structures are still clearly visible in the regime of small regularization parameter values, in particular in the non-sparse parts of the image.

Fig. 4 (a) Screenshot of a series of SIM image reconstructions with generalized Wiener filtering with signal energy (p = 0) regularization and Lukosz bound apodization ( Media 3). Shown is the reconstruction for κ = 5.62 × 10⁻⁴. (b) Screenshot of a series of SIM image reconstructions with generalized Wiener filtering with signal gradient energy (p = 1) regularization and Lukosz bound apodization ( Media 4). Shown is the reconstruction for κ = 10⁻³. (c) Screenshot of a series of SIM image reconstructions with integrated Lukosz filtering with signal energy (p = 0) regularization ( Media 5). Shown is the reconstruction for κ = 10⁻¹. (d) Screenshot of a series of SIM image reconstructions with integrated Lukosz filtering with signal gradient energy (p = 1) regularization and Lukosz bound apodization ( Media 6). Shown is the reconstruction for κ = 1.78 ×10⁻¹. The regularization parameters for the different screenshots have been manually optimized.

Download Full Size | PDF

Additional simulations were performed in order to quantify the sharpness and the noise enhancement effect. The signal-to-noise ratio (SNR) of uniform objects, as well as the full-width at half-maximum (FWHM) and the peak signal-to-noise ratio (PSNR) of point objects, have been computed as a function of the regularization parameter κ, for different values of the side band height parameter (s = 1, s = 2, and s = 3) for both signal energy regularization (p = 0) and signal gradient energy (p = 1) regularization for generalized Wiener filtering with Lukosz bound apodization and for integrated Lukosz filtering. The SNR was measured over the large uniform rectangle in the top left of resolution chart and the PSNR was computed over 64 reconstructed peak signals of point sources for each case.

Fig. 5 shows the resulting graphs, where the SNR and PSNR are normalized to the values obtained for widefield imaging with the same number of collected photons. The FWHM for low values of the regularization parameter is about 0.28λ/NA, which is a factor of about 1.8 smaller than the widefield FWHM-value 0.51λ/NA. The SNR in that regime of regularization parameters is lower than the widefield value, implying that the resolution improvement is obtained at the expense of noise enhancement. The PSNR starts out at a higher value than the widefield case because the narrower spot already implies an increase in peak signal. The FWHM and both the SNR and PSNR increase with increasing regularization parameter, indicating a trade-off between resolution improvement and noise suppression. The exception is the case of integrated Lukosz-filtering with p = 0 regularization, for which violations of the Lukosz bound may occur for side band height parameters larger than approximately 1.6, as discussed in teh appendix. The same trade-off is apparent from the effect of the side band height parameter s. The FWHM and both the SNR and PSNR increase with decreasing s. Signal gradient energy regularization (p = 1) is less sensitive to variations in the side band height parameter and shows flatter curves in the regime of low regularization parameters than signal energy regularization (p = 0), and is therefore the more favourable regularization from the point of view of robustness. Generalized Wiener filtering with Lukosz apodization is to be prefered over integrated Lukosz filtering for p = 0 regularization due to the Lukosz bound violation for large side band height parameters for the integrated Lukosz filtering method.

Fig. 5 Simulated image quality measures as a function of regularization parameter for different side band heights. The columns show the FWHM, the SNR of uniform image patches, and the peak SNR of single spots for generalized Wiener filtering with Lukosz bound apodization with p = 0 (signal energy) regularization ((a)–(c)), for Lukosz apodization with p = 1 (signal gradient energy) regularization ((d)–(f)), integrated Lukosz filtering with p = 0 regularization ((g)–(i)), and for integrated Lukosz filtering with p = 1 regularization ((j)–(l)). The SNR and PSNR are normalized by the values for the corresponding widefield images.

Download Full Size | PDF

The differences between the two filtering approaches for p = 1 regularization are much smaller. A small shift in the regime of regularization parameter values where the resolution-noise-suppression trade-off occurs can be recognized. In addition, integrated Lukosz filtering gives rise to a slightly lower SNR and PSNR curves, due to an effective lower regularization at high spatial frequencies (the regularization term is effectively weighted by the square of the Lukosz bound OTF). This drawback is compensated by a more natural noise spectrum, i.e. the small additional amount of high frequency noise visually masks the low frequency noise structure (see also Fig. 4).

Different apodization functions than the Lukosz bound have been in use so far. A tri-angular apodization function [13], for example, results in a bit sharper image but also in more visible noise artefacts, as it does not comply with the Lukosz bound close to the cut-off. An alternative is the approach based on the distance transform aplied to the support of the OTF [31]. Here the apodization has the form d (v)^ζ, where d (v) is the (normalized) Euclidean distance from the point v to the cutoff of the SIM-OTF and where ζ is an exponent. This approach may be expected to give similar results as the Lukosz bound apodization provided the exponent ζ is a bit higher than one.

The data misfit term of the TM-functional is often interpreted as the logarithm of the likelihood for the reconstructed ground truth image given the measured image data corrupted by noise. A least squares type of data misfit term then corresponds to uniform Gaussian noise. Shot noise, however, is governed by the Poisson distribution and should therefore give rise to a data misfit term in the TM-functional that is different from least squares. It also raises the question if the root cause of the observed noise structure lies in the mismatch between the nature of the noise source and the character of the data misfit term in the TM-functional. To that end we have performed additional simulations where uniform Gaussian noise, e.g. camera readout noise, is added to the raw images. It turns out that the same noise structures are observed as for the shot noise simulations, indicating that the precise nature of the noise source is not particularly relevant for the structured noise artefact.

Finally, we note that the post-filtering technique, in which a filtering operation is applied to the weighted sum of shifted bands reconstruction gives rise to an amount of noise enhancement intermediate between the non-filtered weighted sum of shifted bands reconstruction and the full-blown generalized filtering approach, both for the generalized Wiener-filtering approach followed by Lukosz bound apodization and for the integrated Lukosz bound filtering approach.

6. Experiments

The effects of regularization and apodization have also been investigated using experimental data. We have imaged the chromatin distribution in HDLM-2, which is a Hodgkins lymphoma cell line [33], using a Zeiss Elyra SIM microscope. The cells were fixed on the microscope slide in formaldehyde and the chromatin was counterstained with DAPI (4′,6-diamidino-2-phenylindole). An excitation wavelength of 405 nm was used, the passband of the filtercube was 420–480 nm for the emitted light. A 63X objective lens with numerical aparture of 1.40 was used with immersion oil with refractive index of 1.518. The images were captured on an Andor EM-CCD iXon 885 camera with a pixel size of 8 μm. A tube lens with additional zoom factor 1.6 was used giving a projected pixel size of 79 nm in the object plane and a projected pixel size in the reconstructed SIM images of 43 nm. The projected pixel size of the captured images corresponds to Nyquist sampling for an emission wavelength of 444 nm. The reciprocal pitch of the illumination pattern was calculated from the data using the methods of refs. [31, 34] to be 0.74 times the cut-off frequency of the lens. We have performed additional experiments on actin filaments in bovine pulmonary artery endothelial (BPAE) cells. The F-actin in these cells was stained with Alexa488, the slide was a Molecular Probes prepared slide (F36924). Imaging and reconstruction settings were the same as for the HDLM-2 cells, except the excitation wavelength was 488 nm and the passband of the filtercube was 495–575 nm.

Image reconstruction of a single 2D-slice was done by first applying the phase shift estimation methods developed by Wicker and Heintzmann [31,34] for extracting the spatial frequency bands and then combining and filtering the bands using the Tikhonov-Miller filtering approach with the currently proposed Lukosz bound as apodization function. In the phase shift estimation only the 0th and ±2nd order bands were extracted so as to emulate a true 2D-SIM data acquisition. Fig. 6 and Fig. 7 shows the resulting reconstructions for the two experiments as a function of the regularization parameter for the different filtering and regularization methods also shown in Fig. 4. The other reconstruction parameters were the data misfit weighting parameters α = 0.95 and σ = 1, which is at 50% of the lateral cut-off frequency of the lens, and side band height parameter s = 3. The filters with this Gaussian data misfit term are described in Eq. (24). The relatively high value for the side band height is used to compensate for the loss of high-frequency image content due to the data misfit weighting. The experimental data indicates that Lukosz bound apodization with p = 0 (signal energy) regularization seems to give the best results. The reconstructions appear free from serious artefacts for all values of the regularization parameter, thus providing robust high-resolution images of the structures. Both cases with p = 1 (signal gradient energy) regularization give comparable results, whereas integrated Lukosz filtering with p = 0 (signal energy) regularization gives rise to artefacts for relatively high values of the regularization parameter, in agreement with the simulation results. Fig. 8 and Fig. 9 show a comparison between a widefield image, a weighted sum of shifted bands SIM image and a filtered SIM image. The SIM images shows better detail compared to the widefield image. The weighted sum of bands image, however, has a relatively high level of noise and shows signs of the honeycomb artefact. These flaws are largely suppressed in the filtered SIM image.

Fig. 6 Reconstructions from experimental data of the chromatin distribution in a HDLM-2 cell. (a) Screenshot of a series of reconstructions with generalized Wiener filtering with p = 0 (signal energy) regularization and Lukosz bound apodization ( Media 7). Shown is the reconstruction for κ = 3.16. (b) Screenshot of a series of SIM image reconstructions with generalized Wiener filtering with p = 1 (signal gradient energy) regularization and Lukosz bound apodization ( Media 8). Shown is the reconstruction for κ = 1.78 ×10⁻¹. (c) Screenshot of a series of SIM image reconstructions with integrated Lukosz filtering with p = 0 (signal energy) regularization ( Media 9). Shown is the reconstruction for κ = 3.16 × 10⁻². (d) Screenshot of a SIM image reconstruction with integrated Lukosz filtering with p = 1 (signal gradient energy) regularization and Lukosz bound apodization ( Media 10). Shown is the reconstruction for κ = 5.62 × 10⁻³. The regularization parameters for the different screenshots have been manually optimized. The intensities are in arbitrary units, the scale bars are 5 μm.

Download Full Size | PDF

Fig. 7 Reconstructions from experimental data of actin in BPAE cells. (a) Screenshot of a series of reconstructions with generalized Wiener filtering with p = 0 (signal energy) regularization and Lukosz bound apodization ( Media 11). Shown is the reconstruction for κ = 1 ×10⁻¹. (b) Screenshot of a series of SIM image reconstructions with generalized Wiener filtering with p = 1 (signal gradient energy) regularization and Lukosz bound apodization ( Media 12). Shown is the reconstruction for κ = 1.78 ×10⁻¹. (c) Screenshot of a series of SIM image reconstructions with integrated Lukosz filtering with p = 0 (signal energy) regularization ( Media 13). Shown is the reconstruction for κ = 3.16 ×10⁻². (d) Screenshot of a SIM image reconstruction with integrated Lukosz filtering with p = 1 (signal gradient energy) regularization and Lukosz bound apodization ( Media 14). Shown is the reconstruction for κ = 5.62 ×10⁻¹. The regularization parameters for the different screenshots have been manually optimized. The intensities are in arbitrary units, the scale bars are 5 μm. The cutouts have been linearly stretched for easier visual inspection.

Download Full Size | PDF

Fig. 8 (a) Reconstructed widefield (WF) image of the 2D-slice of the chromatin distribution in a HDLM-2 cell. (b) Power spectrum (on a logarithmic scale) of the reconstructed widefield image. (c) SIM reconstruction with the unfiltered sum of bands method, Eq. (13). (d) Power spectrum (on a logarithmic scale) of the unfiltered sum of bands method image. (e) One of the reconstructed SIM images of Fig. 6 (Lukosz bound apodization with p = 0 signal energy regularization for side band height parameter s = 3 and regularization parameter κ = 3.16) (f) Power spectrum (on a logarithmic scale) of the reconstructed SIM image. The scale of the Fourier space images ((b) and (d, f)) is different due to the difference in number of pixels, consequently the pass-band of the objective lens (the red circles) have a different radius in these figures. The real space images ((a) and (c,e)) are linearly stretched for maximum contrast, the SIM power spectrum (d,f) has been clipped outside the support of the OTF for visibility purposes, scale bars are 5 μm.

Download Full Size | PDF

Fig. 9 (a) Reconstructed widefield (WF) image of the 2D-slice of actin in BPAE cells. (b) Power spectrum (on a logarithmic scale) of the reconstructed widefield image. (c) SIM reconstruction with the unfiltered sum of bands method, Eq. (13). (d) Power spectrum (on a logarithmic scale) of the unfiltered sum of bands method image. (e) One of the reconstructed SIM images of Fig. 7 (Lukosz bound apodization with p = 0 signal energy regularization for side band height parameter s = 3 and regularization parameter κ = 1 ×10⁻¹) (f) Power spectrum (on a logarithmic scale) of the reconstructed SIM image. The scale of the Fourier space images ((b) and (d, f)) is different due to the difference in number of pixels, consequently the pass-band of the objective lens (the red circles) have a different radius in these figures. The real space images ((a) and (c,e)) are linearly stretched for maximum contrast, the cutouts are linearly contrast stretched in the insets as well. The SIM power spectrum (d,f) has been clipped outside the support of the OTF for visibility purposes, scale bars are 5 μm.

Download Full Size | PDF

7. Discussion

In summary, we have analyzed regularization and apodization in image reconstruction for 2D-SIM based on Tikhonov-Miller filtering and using the Lukosz bound as apodization function. The unfiltered image reconstruction method based on the weighted addition of spatial frequency bands gives rise to both resolution improvement and higher noise levels for increasing values of the relative side band height parameter. The optimum of this trade-off is found for relative side band heights in the range 1–2. The noise enhancement can be largely suppressed by the introduction of Tikhonov-Miller filters. In particular, the generalized filtering technique allows for the best noise suppression. The Lukosz bound can be used as a post-processing apodization filter or can be integrated into the Tikohonov-Miller framework. The latter option results in a slight increase in the high frequency range of the residual noise structure. Practical 2D-reconstructions of thick objects generally suffer from out-of-focus background. The honeycomb artefacts, which result from this out-of-focus light, can be suppressed by a data misfit weight function. Reasonable results can be obtained with these filters in combination with elevated side band height parameters (in the range 3–6). No clipping of negative pixel values to zero was needed for any of the reconstruction methods based on the Lukosz bound, a practice that is a default setting in some commercial systems.

There are a number of ways to automatically set the regularization parameter [35]. Generalized Cross-Validation (GCV) is often proposed as a suitable method for finding the ‘optimum’ regularization parameter in the field of image restoration [35, 36]. According to this method the best value for the regularization parameter is the one that minimizs the cross-validation (CV) function. This function is defined as the scaled version of the data misfit term in the TM-functional. The CV-function for generalized filtering is:

C V_{gen} (κ) = \frac{\int \sum_{lm} {| {\hat{B}}_{lm} (v) {\hat{I}}_{gen} (v) - s_{m} {\hat{I}}_{band, lm} (v + R_{l} \cdot q_{m}) |}^{2} d v}{{[\int (\sum_{lm} {\hat{B}}_{lm} (v) {\hat{F}}_{lm} (v) - 1) d v]}^{2}} .

Note that the CV-function is slightly different when the data misfit weight function is used. We have noted that the regularization parameter κ_GCV following from the generalized cross-validation (GCV) method does not always correspond to what intuitively would be the ‘best’ image. This way of automatically setting the regularization paramater must therefore be exercised with care.

We envision that the current research can be extended in a number of ways. First, an extension to a full 3D reconstruction scheme seems logical to explore. A second extension warranted by the results presented here is a systematic study into noise propagation in SIM image reconstruction. Such a study should incorporate effects of shot noise, readout noise and phase stepping jitter. Finally, different ways of regularization, better suited to suppressing typical SIM artefacts, may be explored.

Appendix

This integrated Lukosz filtering approach gives rise to an effective OTF that satisfies the Lukosz bound for all values of the regularization parameter, with the exception of the signal energy regularization case (p = 0). The normalization of the OTF (must be equal to one for v = 0) implies that the r.h.s. of Eq. (30) must be multiplied with a normalization factor:

η_{norm} = 1 + \frac{κ}{\sum_{lm} s_{m}^{2} w_{m}^{2} {| \hat{H} (R_{l} \cdot q_{m}) |}^{2}} .

Now the effective OTF is no longer necessarily a decreasing function of the regularization parameter κ. It follows that violations of the Lukosz bound can then occur, as the effective OTF is equal to the Lukosz bound for κ = 0. The requirement that the effective OTF must satisfy the Lukosz bound leads to the constraint that the function:

z (v) = \frac{\sum_{lm} s_{m}^{2} w_{m}^{2} {| \hat{H} (v + R_{l} \cdot q_{m}) |}^{2}}{\hat{Λ} {(v)}^{2}} .

must be smaller than its value for v = 0 for all v. It may be expected that, whenever this constraint is not satisfied, this occurs for spatial frequencies at the side band peaks, as the function z (v) has local maxima there. For these spatial frequencies it may be deduced that the constraint is satisfied provided the side band height parameter is sufficiently small:

s \leq s_{max} = \frac{1}{w} \sqrt{\frac{3 Λ_{s}^{2} - H_{s}^{2}}{1 - 2 (3 Λ_{s}^{2} - 1) H_{s}^{2}}},

with H_s and Λ_s equal to the widefield OTF and to the Lukosz bound at the side-peak spatial frequency, respectively. For a side band peak at 90% of the wide field cut-off (H_s = 0.037 and Λ_s = 0.45) and a maximum modulation depth of the illumination (w = 1/2) we find s_max = 1.56.

Acknowledgments

Rainer Heintzmann and Kai Wicker are thanked for making their extensive SIM processing software package available to us, for advice and for suggesting the use of a digital resolution chart in simulations. Adriaan Houtsmuller is acknowledged for support in acquiring experimental data.

References and links

1. S. W. Hell, “Microscopy and its focal switch,” Nat. Methods 6, 24–32 (2009). [CrossRef] [PubMed]

2. L. Schermelleh, R. Heintzmann, and H. Leonhardt, “A guide to super-resolution fluorescence microscopy,” J. Cell Biol. 190, 165–175 (2012). [CrossRef]

3. W. Lukosz, “Optical systems with resolving powers exceeding the classical limit,” J. Opt. Soc. Am. 56, 1463–1471 (1966). [CrossRef]

4. M. A. A. Neil, R. Juskaitis, and T. Wilson, “Method of obtaining optical sectioning by using structured light in a conventional microscope,” Opt. Lett. 22, 1905–1907 (1997). [CrossRef]

5. M. A. A. Neil, A. Squire, R. Juskaitis, P. I. H. Bastiaens, and T. Wilson, “Wide-field optically sectioning fluorescence microscopy with laser illumination,” J. Microsc. 197, 1–4 (2000). [CrossRef] [PubMed]

6. R. Heintzmann and C. Cremer, “Laterally modulated excitation microscopy: improvement of resolution by using a diffraction grating,” Proc. SPIE 3568, 185–196 (1999). [CrossRef]

7. G. E. Cragg and P. T. C. So, “Lateral resolution enhancement with standing evanescent waves,” Opt. Lett. 25, 46–48 (2000). [CrossRef]

8. M. G. L. Gustafsson, “Surpassing the lateral resolution limit by a factor of two using structured illumination microscopy,” J. Microsc. 198, 82–87 (2000). [CrossRef] [PubMed]

9. J. T. Frohn, H. F. Knapp, and A. Stemmer, “True optical resolution beyond the Rayleigh limit achieved by standing wave illumination,” Proc. Natl. Acad. Sci. U.S.A. 97, 7232–7236 (2000). [CrossRef] [PubMed]

10. R. Heintzmann, T. Jovin, and C. Cremer, “Saturated patterned excitation microscopy - a concept for optical resolution improvement,” J. Opt. Soc. Am. B 19, 1599–1609 (2002). [CrossRef]

11. M. G. L. Gustafsson, “Nonlinear structured-illumination microscopy: Wide-field fluorescence imaging with theoretically unlimited resolution,” Proc. Natl. Acad. Sci. U.S.A. 102, 13081–13086 (2005). [CrossRef] [PubMed]

12. L. Wang, M. C. Pitter, and M. G. Somekh, “Wide-field high-resolution structured illumination solid immersion fluorescence microscopy,” Opt. Lett. 36, 2794–2796 (2011). [CrossRef] [PubMed]

13. M. G. L. Gustafsson, L. Shao, P. M. Carlton, C. J. R. Wang, I. N. Golubovskaya, W. Z. Cande, D. A. Agard, and J. W. Sedat, “Three-dimensional resolution doubling in wide-field fluorescence microscopy by structured illumination,” Biophys. J. 94, 4957–4970 (2008). [CrossRef] [PubMed]

14. P. Kner, B. B. Chhun, E. R. Griffis, L. Winoto, and M. G. L. Gustafsson, “Super-resolution video microscopy of live cells by structured illumination,” Nat. Methods 6, 339–342 (2009). [CrossRef] [PubMed]

15. L. Shao, P. Kner, E. H. Rego, and M. G. L. Gustafsson, “Super-resolution 3d microscopy of live whole cells using structured illumination,” Nat. Methods 12, 1044–1046 (2011). [CrossRef]

16. R. Fiolka, M. Beck, and A. Stemmer, “Structured illumination in total internal reflection fluorescence microscopy using a spatial light modulator,” Opt. Lett. 33, 1629–1631 (2008). [CrossRef] [PubMed]

17. O. Mandula, M. Kielhorn, K. Wicker, G. Krampert, I. Kleppe, and R. Heintzmann, “Line scan - structured illumination microscopy super-resolution imaging in thick fluorescent samples,” Opt. Express 20, 24167–24174 (2012). [CrossRef] [PubMed]

18. G. M. P. van Kempen, L. J. van Vliet, P. Verveer, and H. van der Voort, “A quantitative comparison of image restoration methods for confocal microscopy,” J. Microsc. 185, 354–365 (1997). [CrossRef]

19. P. Pankajakshan, B. Zhang, L. Blanc-Feraud, Z. Kam, J. C. Olivo-Marin, and J. Zerubia, “Parametric blind deconvolution for confocal laser scanning microscopy,” Conf. Proc. IEEE Eng. Med. Biol. Soc. 2007, 6532–6535 (2007).

20. G. M. P. van Kempen and L. J. van Vliet, “Background estimation in nonlinear image restoration,” J. Opt. Soc. Am. A 17, 425–434 (2000). [CrossRef]

21. N. Dey, L. Blanc-Feraud, C. Zimmer, P. Roux, Z. Kam, J. C. Olivo-Marin, and J. Zerubia, “Richardson-Lucy algorithm with total variation regularization for 3D confocal microscope deconvolution,” Microsc. Res. Tech. 69, 260–266 (2006). [CrossRef] [PubMed]

22. C. Berenstein and E. Patrick, “Exact deconvolution for multiple convolution operators–an overview, plus performance characterizations for imaging sensors,” Proc. IEEE 78, 723–734 (1990). [CrossRef]

23. L. P. Yaroslavsky and H. J. Caulfield, “Deconvolution of multiple images of the same object,” Appl. Opt. 33, 2157–2162 (1994). [CrossRef] [PubMed]

24. S. A. Shroff, J. R. Fienup, and D. R. Williams, “Phase-shift estimation in sinusoidally illuminated images for lateral superresolution,” J. Opt. Soc. Am. A 26, 413–424 (2009). [CrossRef]

25. W. Lukosz, “Übertragung nicht-negativer signale durch lineare filter,” J. Mod. Opt. 9, 335–364 (1962).

26. W. Lukosz, “Properties of linear low-pass filters for nonnegative signals,” J. Opt. Soc. Am. 52, 827–829 (1962). [CrossRef]

27. P. T. C. So, H.-S. Kwon, and C. Y. Dong, “Resolution enhancement in standing-wave total internal reflection microscopy: a point-spread-function engineering approach,” J. Opt. Soc. Am. A 25, 1319–1329 (2008).

28. M. G. Somekh, K. Hsu, and M. C. Pitter, “Resolution in structured illumination microscopy: a probabilistic approach,” J. Opt. Soc. Am. A 25, 1319–1329 (2008). [CrossRef]

29. A. N. Tikhonov and V. A. Arsenin, Solution of Ill-posed Problems (Winston-Wiley, 1977).

30. N. Wiener, Extrapolation, Interpolation, and Smoothing of Stationary Time Series (Wiley, 1949).

31. K. Wicker, O. Mandula, G. Best, R. Fiolka, and R. Heintzmann, “Phase optimization for structured illumination microscopy,” Opt. Express 21, 2032–2049 (2013). [CrossRef] [PubMed]

32. C. Luengo Hendriks, B. Rieger, M. van Ginkel, G. M. P. van Kempen, and L. J. van Vliet, “DIP-image: a scientific image processing toolbox for MATLAB,” (1999-). Delft University of Technology, http://www.diplib.org/.

33. H. G. Drexler, G. Gaedicke, M. S. Lok, V. Diehl, and J. Minowada, “Hodgkin’s disease derived cell lines HDLM-2 and L-428: comparison of morphology, immunological and isoenzyme profiles,” Leuk. Res. 10, 487–500 (1986). [CrossRef]

34. K. Wicker, PhD thesis, King’s College London, (2010).

35. G. M. P. van Kempen and L. J. van Vliet, “The influence of the regularization parameter and the first estimate on the performance of Tikhonov regularized non-linear image restoration algorithms,” J. Microsc. 198, 63–75 (2000). [CrossRef] [PubMed]

36. N. P. Galatsanos and A. K. Katsaggelos, “Methods for choosing the regularization parameter and estimating the noise variance in image restoration and their relation,” IEEE Trans. Image Process. 1, 322–336 (1992). [CrossRef] [PubMed]

Image filtering in structured illumination microscopy using the Lukosz bound

Abstract

1. Introduction

2. The Lukosz bound

3. Image formation in SIM

3.1. Image acquisition and phase shifting

3.2. Weighted sum of shifted bands reconstruction

3.3. Filtering approaches

3.4. Data misfit weight function

4. The Lukosz bound in SIM filtering

4.1. 2D SIM Lukosz bound

4.2. Lukosz bound filtering

5. Simulation results

6. Experiments

7. Discussion

Appendix

Acknowledgments

References and links

Supplementary Material (14)

Cited By

Figures (9)

Equations (34)

Optics Express