Sted microscope method for improving signal-to-noise ratio under low-photon counting imaging conditions

Abstract

The invention relates to a method for generating a result image. The method comprises the steps of: acquiring a STED image of a sample, the STED image comprising pixels; calculating Fourier coefficients of the time of arrival for the pixels of the image, resulting in real part coefficients representing a first image and imaginary part coefficients representing a second image; deriving an intensity image from the STED image; applying a spatial filter to the first image, the second image and the intensity image, resulting in respective filtered images; calculating an image G based on the filtered first image and the filtered intensity image; calculating an image S based on the filtered second image and the filtered intensity image; calculating the result image from the two images G and S.

Description

Background Technology

[0001] In STED (stimulated emission depletion) fluorescence microscopy, as described, for example, in EP 0801 759 B1, at least one image of a sample stained with a fluorescent molecule or dye is obtained by exciting the fluorescent molecule or dye with a light pulse in the focal region of the microscope objective. Furthermore, de-excitation light is applied to the sample to de-excite the region outside the center of the focal region, or the focal region itself. Therefore, this type of microscopy can produce images with spatial resolution below the usual optical diffraction limit. De-excitation can originate from a pulsed light source or a continuous wave (CW) light source. Both of these sources are typically lasers, referred to as excitation and depletion lasers.

[0002] US 9,551,658 B2 describes a method and apparatus for extracting light pulses derived from a sample with better spatial resolution using time-gated detection. De-excitation light leads to a shortened fluorescence lifetime. Therefore, the time decay of the fluorescence signal outside the center of the focal region is shorter than the time decay at the center. In this case, the lifetime τ is... Given, where k nr It is the nonradiative deexcitation rate, k r It is the radiative deexcitation rate, k STED The additional nonradiative de-excitation rate is introduced by the depletion of the laser. The time gate excludes the signal within a short period after the excitation pulse. The remaining signal is collected. In this time-gated method, the collection of light originating from the center of the focal region is also incomplete. Signal loss occurs from photons arriving during the gate-off period, resulting in a lower signal-to-noise ratio than the lower-resolution non-gated image.

[0003] Publication: Lanzanò et al. (2015), Encoding and decoding spatio-temporal information for super-resolution microscopy, Nature Communications, 6:6701, 10.1038 / ncomms7701, describes a SPLIT method that preserves the signal of early photons by applying a phasor method. Here, the lifetime decay of each pixel is converted into a polar coordinate map, called a phasor map. In this map, the coordinates of each pixel, called a phasor, are a linear combination of the coordinates of three components: the high-resolution signal, the low-resolution signal, and the background. The weights of the high-resolution components are used as the intensity of the resulting image and contain information from all detected photons. Although the signal of all photons is used, this method has the limitation of requiring a large number of photons to be acquired. Under low photon count conditions, the distribution in the phasor map is so wide that it is impossible to achieve an improvement in the signal-to-noise ratio.

[0004] To achieve a narrower distribution, phasor coordinates must be filtered. A chapter in this book, Digman et al. (2012), *The phasor approach to fluorescence lifetime imaging: exploiting phasor linear properties, in Fluorescence Lifetime Spectroscopy and Imaging: Principles and Applications in Biomedical Diagnostics* (edited by L. Marcu, PMW French, and DS Elson), CRC Press, describes a method for filtering phasor maps with the aim of preserving spatial information through the use of a median filter. However, in the case of super-resolution, even this filtering method can lead to an unnecessarily decreased spatial resolution.

[0005] In digital image processing, more advanced filtering methods are known. For example, the publication Sendur, L. & Selesnick, I. (2003), Bivariate Shrinkage with Local Variance Estimation, Signal Processing Letters, IEEE 9,438–441, 10.1109 / LSP.2002.806054, describes an advanced denoising filter based on wavelet transform, which better preserves spatial resolution when applied to images. Attempts to apply this filter or similar filters to the phasor components in the same manner as median filters have failed due to the strong variations in noise in the phasor component images. Summary of the Invention

[0006] The object of this invention is to provide a method for increasing the signal-to-noise ratio in STED images under low photon count conditions while maintaining the spatial resolution of the gated STED method. This object is achieved by a method comprising at least the following steps:

[0007] Obtain the STED image of the sample, which includes pixels;

[0008] Calculate the Fourier coefficients or photon arrival time data for each pixel in the image to obtain the real and imaginary coefficients.

[0009] The real part coefficients represent the first image, and the imaginary part coefficients represent the second image;

[0010] Intensity images are derived from STED images;

[0011] Spatial filters are applied to the first image, the second image, and the intensity image to obtain their respective filtered images;

[0012] Image G is calculated based on the filtered first image and the filtered intensity image;

[0013] Image S is calculated based on the filtered second image and the filtered intensity image; and

[0014] The resulting image is calculated based on two images, G and S.

[0015] Preferably, acquiring an STED image includes the following steps:

[0016] The sample is irradiated with pulsed excitation light, wherein the pulsed excitation light includes a focal point when irradiating the sample;

[0017] The sample is illuminated with excitation light surrounding the central region of the focal point.

[0018] Typically, two-dimensional (and sometimes three-dimensional) STED images of a sample or portion of a sample are acquired using a scanning microscope with STED functionality. For example, an acquired two-dimensional STED image comprises pixels in the x and y directions (e.g., 512 x 512 pixels), depending on the acquisition parameters of the scanning microscope. For each pixel, the photon arrival time relative to the laser pulse is recorded. The fluorescence lifetime of typical fluorescent dyes ranges from 1 to 4 ns. Many STED dyes have multiple lifetime compositions. There is no or weak deexcitation light at the center of the focal region, and the natural lifetime is observed there. Outside the center of the focal region, the lifetime decreases / is lower, but typically below 1 ns, for example, in areas where deexcitation light is applied, due to the deexcitation process. The decrease in lifetime outside the center of the focal region is inversely proportional to the depletion energy from the depletion laser. For a typical STED image, approximately 100 photon counts are recorded in the bright pixels. The Fourier coefficients of the photon arrival time are calculated pixel-by-pixel, preferably for all pixels of the image. Because each pixel has real and imaginary Fourier coefficients, the real (Fourier) coefficients can be considered as the first image, and the imaginary (Fourier) coefficients as the second image. However, it is not important whether the first or second image actually forms, for example, a two-dimensional image. The first and second images can be represented or considered as data structures with, for example, only linear or other suitable data storage arrangements. A similar situation applies to at least one of the following: the resulting image, the filtered image, the intensity image, image G, and image S.

[0019] In a preferred embodiment, acquiring the STED image further includes registering the arrival time of the detected fluorescent photons relative to the excitation pulse for each pixel.

[0020] Alternatively or additionally, the arrival time is registered through at least two time gates.

[0021] Preferably, the calculation of Fourier coefficients includes the calculation of the nth Fourier coefficient, where n is a pre-selected value (a natural number).

[0022] Preferably, the intensity image is obtained from an STED image or calculated from intensity values ​​derived from the number of photons registered for a pixel.

[0023] In a preferred embodiment, image G is calculated based on dividing the filtered first image by the filtered intensity image.

[0024] In another preferred embodiment, image S is calculated based on dividing the filtered second image by the filtered intensity image.

[0025] In a further preferred embodiment, the spatial filter is at least one of anisotropic diffusion, total variation denoising, maximum entropy denoising, nonlocal mean, and small wavelet denoising.

[0026] The above objectives are further achieved by means of a device suitable for performing the above methods.

[0027] Furthermore, the aforementioned objective is also achieved by a computer program that executes the above-described method. The first and / or second images in the sense of this application can be represented in the form of a two-dimensional array of corresponding Fourier coefficients.

[0028] In a particularly preferred embodiment of the invention, a filter / spatial filter of the first image, the second image, and the intensity image is applied before dividing the Fourier coefficients by the filtered intensity. Previously, the filter was applied after the method step of performing division on the image. The phasor coordinates of adjacent pixels (after division) have similar values, while the original Fourier coefficients vary depending on the intensity, especially under low photon count conditions. In this context, the term phasor coordinates is used in a manner similar to that described by Lanzanò et al. For a smaller number of photon counts, the Fourier coefficients of adjacent pixels show greater variation than for a higher photon count. The different variations in similar average values ​​contradict the prerequisites of advanced filtering methods. In the prior art, a median filter is chosen, which is reasonable when the photon counts in adjacent pixels show little variation. However, the STED image shows a high degree of variation that must be preserved. According to the invention, when the image is filtered before performing division on the image, prerequisites for advanced filtering are given. Furthermore, when performing division on a filtered image, there is a bias that the filtered image may contain information from different regions with different weights, and therefore the division step may result in a more widespread distribution of phasor coordinates.

[0029] However, the results obtained through the process according to the particularly preferred embodiment are unexpectedly superior to classical / previous methods, to the point that the objective is practically achievable. It has been found that the advanced filtering algorithm is more sensitive to strong changes in statistical data among neighboring pixels. Attached Figure Description

[0030] Figure 1 illustrates the effect of STED de-excitation on lifetime decay, and the influence of... Figure 1a The central area of ​​the focus, Figure 1b A typical lifetime decay phasor diagram of data in the central outer region of the image. Figure 1c The typical choices for lifetime decay and time gates for the combination are shown. The positions of these components in the phasor diagram are shown in Figure 1.d.

[0031] Figure 2 shows the phasor distribution of twenty image pixels when different filtering methods are applied according to a preferred embodiment of the present invention: Figure 2a No filter Figure 2b It is a 5x5 median filter. Figure 2c It is a 7x7 median filter, and Figure 2d It's a wavelet filter. The data used is simulated based on the average of ten photon counts per pixel.

[0032] Figure 3 shows the target Figure 3b The best filter in Figure 3c There is no filter in it. Figure 3d The median filter 7x7 and Figure 3e The filter effect of the filter according to a preferred embodiment of the present invention with respect to undesirable smooth edges. Detailed Implementation

[0033] In addition to the standard STED microscope typically implemented in a confocal laser scanning microscope, a device for detecting the photon arrival time relative to the excitation pulse of the excitation light is mandatory for the method according to the invention. This device may, for example, comprise a recording apparatus based on the time-dependent single-photon counting principle (TCSPC principle) or a time-gating recording apparatus.

[0034] The recorded data can be represented as the number of photons N detected at arrival time t for pixel i: Ni(t). For each pixel, the Fourier coefficients Freal and Fimag are calculated:

[0035]

[0036] The integer value n is the pre-selected order of the Fourier coefficients. Typically, the first order (1) is chosen, i.e., n = 1. f is the pulse repetition rate of the excitation light. Furthermore, for each pixel, the intensity I is calculated:

[0037]

[0038] In the next step, spatial filters, particularly denoising filters, are used to filter the three images Freal, Fimag, and I. Preferably, several advanced denoising algorithms that preserve spatial resolution can be used in this step. For example, algorithms such as anisotropic diffusion, total variation denoising, maximum entropy denoising, and non-local-means and wavelet denoising can be applied.

[0039] In a preferred embodiment, a wavelet filter is selected. The basic principle is similar to that of a Fourier filter used for denoising. Gaussian noise causes an additional shift in the amplitude of all Fourier coefficients. The image is transformed to Fourier space. The shift is subtracted, and the image is back-transformed. This process is also called shrinkage. The difference between wavelet space and Fourier space is that wavelet coefficients represent frequency information for regions in the image. For wavelet denoising, the same principle of coefficient shrinkage is used. The simplest method is to shrink all coefficients by the same shift. More advanced methods that can be applied to this invention have been developed, for example: Sendur, L. & Selesnick, I., (2003), Bivariate Shrinkage with Local Variance Estimation, Signal Processing Letters, IEEE. Vol. 9. 438–441, 10.1109 / LSP. 2002. 806054.

[0040] The contraction method assumes a Gaussian distribution for the noise. For the intensity image I, the noise distribution is known to be a Poisson distribution. One method of using a filter designed to convert Poisson-distributed noise to Gaussian-distributed noise is, for example, by transforming the data through the Anscombe transform:

[0041]

[0042] The filter is then applied, and the data is inversely transformed. The standard deviation of the noise after the Anscombe transform is 1.0.

[0043] This situation differs for the two images, Freal and Fimag. After integration, the noise is not perfectly Poisson distributed, but depends on the fluorescence lifetime.

[0044] However, this situation is not as severe as attempts to filter the phasor components. In the first order, it can be assumed that the distribution is almost Poisson, and the noise variance for these images can be estimated using statistically robust methods such as the median absolute deviation (MAD):

[0045] σ=meduan(|X i -median(X i )|) / 0.67449 (4)

[0046] Where σ is the standard deviation, and the median operation is applied to the intensity X of all pixels i. There is a dependency, specifically a linear dependency, between σ and the optimal shrinkage offset. The factor depends on the implementation details of the wavelet transform and can be derived for the specific shrinkage algorithm used. For example, in Sendur, L. & Selesnick, I., the dependency is given in Equation 5, and Part II contains the derivation.

[0047] In a preferred embodiment, all three images, Freal, Fimag, and I, are transformed using the Anscombe transform (3), and then noise is estimated for all three images using the median absolute deviation method (4). The Anscombe-transformed images are transformed to wavelet space. The coefficients are shrunk using values ​​derived from the estimated standard deviation of the noise. The shrunk images are then inversely transformed in wavelet space and subjected to an inverse Anscombe transform. The result is the filtered images Freal_filtered, Fimag_filtered, and I_filtered.

[0048] In the next step, the phasor coefficients for all image pixels are calculated:

[0049] G i =Freal_filtered i / I_filtered i ,S i =Fimag_filtered i / I_filtered i (5)

[0050] In the final step, the resulting image is calculated based on the phasor coefficients. This calculation can be performed using a method similar to the SPLIT method described in Lanzano et al. (2015), Encoding and decoding spatio-temporal information for super-resolution microscopy, Nature Communications, 10.1038 / ncomms7701.

[0051] The following explains why the method according to the present invention outperforms existing techniques for filtering the phasor G and S components. Figure 1a The diagram shows the attenuation of intensity I along the vertical axis 1 and arrival time t along the horizontal axis 2, where Figure 3 represents photons originating from the focal center region of the excitation light. Figure 1bChart 6 shows photons originating from the outer region of the focal point, in which excitation light is applied. The average arrival time of photons originating from the outer region of the focal point is shorter than the average arrival time of photons originating from the focal point. Figure 1c Chart 7 shows all detected photons. In the time-gated STED, the arrival time interval 8 is selected (e.g., by gating) where photons from the central region dominate. Photons detected before this interval are ignored. Therefore, signal 4 (see [reference]) contains information about this central region. Figure 1a ) is missing. Only part 5 of the signal (see Figure 1a This method is used to generate images. Compared to the classic STED method, the temporal gating method produces images with better spatial resolution but higher noise levels.

[0052] The pixel with attenuation originating solely from the central region is located, for example, at position 9 in the phasor diagram (see...). Figure 1d Photons from the outside are located at 10 different positions in the phasor diagram. The idea behind the SPLIT method (prior art) is to apply weights to the intensity of a pixel based on its position in the phasor diagram. This achieves spatial resolution comparable to the time-gated STED method, but with less noise because no photons are lost to be counted.

[0053] Figure 2a Simulated data for the average of 10 photon counts in an image pixel is visualized. The simulation assumes a mixture of 50% photons from the central region with a lifetime of 4n (4 nanoseconds) 11 and 50% photons from the outside with a lifetime of 200 ps (200 picoseconds) 12. Phasor coordinates 13 for 20 pixels are plotted. For high photon counts, the phasor position must be halfway between the positions of the two components 11 and 12, i.e., at 14. For the case of low photon counts without filtering, represented by 15, it is assumed that due to the widespread distribution without filtering, there may be no improvement compared to the gating method.

[0054] Figure 2b The same simulation data is shown when a 5x5 median filter, represented by 16, is applied according to existing technology. For Figure 2c A 7x7 median filter, represented by 17, was used. Finally, Figure 2d The results (denoted by 18) when using a wavelet filter according to an embodiment of the present invention are shown. A 7x7 median filter is required to achieve a comparable signal-to-noise ratio.

[0055] The filter must additionally preserve spatial resolution to be useful in STED applications. A simulation with two mixtures illustrates the filter's behavior along the edge. According to... Figure 3aIn the simulation, in the left half 20, image 19 is created by a mixture of 2 / 3 low-resolution components and 1 / 3 high-resolution components. In the right half 21, 1 / 3 of the low-resolution components and 2 / 3 of the high-resolution components are simulated. The 7 rows of 8 pixels in the central region 22 have been vertically accumulated and weighted with an intensity gradient along a line in the phasor between 4ns and 200ps.

[0056] Figure 3b The optimal conditions for very high photon counts are shown, with an intensity of 23 for an average of 8 pixels. Figure 3c For cases without filtering, in Figure 3d The figures for the 7x7 median filter and the wavelet filter according to the invention in Figure 2e show the intensity for 8 average pixels with 10 counts per pixel.

[0057] The wavelet filter according to the present invention is superior to the median filter according to the prior art. Figure 3d It better maintains spatial resolution. Figure 3e ), and has a considerable signal-to-noise ratio ( Figure 2c and Figure 2d ).

[0058] An advanced filter has been designed to address the similar noise characteristics in adjacent pixels. However, this filter fails when directly applied to phasor components using existing techniques. This is because the noise is significantly altered due to the intensity division in the phasor formula.

[0059] Some or all of the method steps may be performed by (or using) hardware devices, such as processors, microprocessors, programmable computers, or electronic circuits. In some embodiments, one or more of the most important method steps may be performed by such devices.

[0060] Depending on certain implementation requirements, embodiments of the present invention can be implemented in hardware or software. This implementation can be performed using non-transitory storage media, such as digital storage media like floppy disks, DVDs, Blu-ray discs, CDs, ROMs, PROMs, EPROMs, EEPROMs, or FLASH memories, storing electronically readable control signals that cooperate (or are capable of cooperating with) a programmable computer system to perform the corresponding methods. Therefore, the digital storage medium can be computer-readable.

[0061] Some embodiments of the invention include a data carrier having electronically readable control signals, which is capable of cooperating with a programmable computer system to perform one of the methods described herein.

[0062] Typically, embodiments of the present invention can be implemented as a computer program product having program code that, when run on a computer, is operable to perform one of the methods. For example, the program code may be stored on a machine-readable medium.

[0063] Other embodiments include a computer program stored on a machine-readable medium for performing one of the methods described herein.

[0064] In other words, therefore, embodiments of the present invention are computer programs having program code for performing one of the methods described herein when the computer program is run on a computer.

[0065] Therefore, a further embodiment of the invention is a storage medium (or data carrier, or computer-readable medium) including a computer program stored thereon for performing one of the methods described herein when executed by a processor. Data carriers, digital storage media, or recording media are typically tangible and / or non-transient. A further embodiment of the invention is an apparatus as described herein, including a processor and a storage medium.

[0066] Therefore, a further embodiment of the invention represents a data stream or signal sequence for performing one of the methods described herein. The data stream or signal sequence may, for example, be configured to be transmitted via a data communication connection, such as via the Internet.

[0067] Further embodiments include processing means, such as a computer or programmable logic device, configured or adapted to perform one of the methods described herein.

[0068] Further embodiments include a computer on which a computer program for performing one of the methods described herein is installed.

[0069] Further embodiments of the invention include an apparatus or system configured to transmit (e.g., electronically or optically) a computer program for performing one of the methods described herein to a receiver. For example, the receiver may be a computer, mobile device, memory device, etc. For example, the apparatus or system may include a file server for transmitting the computer program to the receiver.

[0070] In some embodiments, a programmable logic device (e.g., a field-programmable gate array) may be used to perform some or all of the functions of the methods described herein. In some embodiments, the field-programmable gate array may cooperate with a microprocessor to perform one of the methods described herein. Generally, these methods are preferably performed by any hardware device.

Claims

1. A method for generating a result image, the method comprising the steps of: a. acquiring a STED image of a sample, the STED image comprising pixels; b. calculating Fourier coefficients of the time of arrival of photons for the pixels of the STED image, resulting in real coefficients representing a first image and imaginary coefficients representing a second image; c. deriving an intensity image from the STED image; d. applying a spatial filter to the first image, the second image and the intensity image, resulting in respective filtered images; e. calculating an image G based on the filtered first image and the filtered intensity image by dividing the filtered first image by the filtered intensity image; f. calculating an image S based on the filtered second image and the filtered intensity image by dividing the filtered second image by the filtered intensity image; and g. calculating the result image based on the two images G and S.

2. The method according to claim 1, wherein acquiring the STED image comprises the steps of: illuminating the sample with a pulsed excitation light having a focal point; and illuminating the sample with a de-excitation light surrounding a central region of the focal point.

3. The method according to claim 1 or 2, wherein acquiring the STED image further comprises: registering the time of arrival of the detected fluorescent photons with respect to the excitation pulse for each pixel separately, and / or registering the time of arrival by at least two time gates.

4. The method according to claim 1 or 2, wherein the calculation of the Fourier coefficients comprises the calculation of the n-th order Fourier coefficients, wherein n is a preselected value. The intensity image is derived from the STED image or is derived from the calculation of intensity values from the number of photons registered for the pixels.

6. The method according to claim 1 or 2, further comprising the steps of:

5. The method of claim 1 or 2, wherein, transforming the first image, the second image and the intensity image having a first noise distribution into transformed images having a second noise distribution, and performing the method steps d-g based on the transformed images.

7. The method according to claim 6, wherein the first noise distribution is a Poisson noise distribution and the second noise distribution is a Gaussian noise distribution.

8. The method according to claim 6, further comprising the step of inverse transforming the filtered first image, the filtered second image and the filtered intensity image after performing step d. The spatial filter is at least one of anisotropic diffusion, total variation denoising, maximum entropy denoising, non-local means and wavelet denoising.

10. A device adapted to perform the method according to any of the preceding claims, wherein the device is a STED microscope or a STED scanning microscope.

9. The method of claim 1 or 2, wherein, 11. A computer program product comprising program code for performing the method according to any of claims 1 to 9 when executed on a processor. ​ ​

Images