Bleach correction ImageJ plugin for compensating the photobleaching of time-lapse sequences

During the capturing of the time-lapse sequence of fluorescently labeled samples, fluorescence intensity exhibits decays. This phenomenon is known as 'photobleaching' and is a widely known problem in imaging in life sciences. The photobleaching can be attenuated by tuning the imaging set-up, but when such adjustments only partially work, the image sequence can be corrected for the loss of intensity in order to precisely segment the target structure or to quantify true intensity dynamics. We implemented an ImageJ plugin that allows the user to compensate for the photobleaching to estimate the non-bleaching condition with choice of three different algorithms: simple ratio, exponential fitting, and histogram matching methods. The histogram matching method is a novel algorithm for photobleaching correction. This article presents details and characteristics of each algorithm based on application to actual image sequences.


Introduction
In biological fluorescence microscopy, cells are irradiated with excitation light that causes the emission of fluorescence from protein markers. This irradiation is necessary to detect the position of proteins but it also has a side effect. Emitted fluorescence gradually decreases by time. This is because, with a certain probability, fluorophores irreversibly lose the ability to fluoresce. This effect is called "photobleaching" and is a widely-known problem in bioimage analysis. Photobleaching not only degrades the visual quality of the results but also interferes with the measurement of molecular kinetics and the quality of the segmentation of target objects. Since bleaching attenuates total signal intensity even when the density of labeled protein is unchanged, precise estimation of the amount of protein or the boundary of the biological structure becomes difficult.
To overcome the problem of photobleaching, improvements can be made either before or after the experiment. Before the experiment, one can tune instruments and sample environment, such as careful choice of the fluorophore, use of anti-fading reagents, decreasing the power of laser irradiation, increasing the detector gain, and increasing the time interval of capturing 1 . After the experiment, time-lapse sequence image data can be processed to compensate for the loss of intensity. To do this, the amount of the fluorescence loss is estimated and then the none-bleached condition is restored by image processing. We call such restoration procedure "bleach correction". Several different algorithms have been developed and used by various researchers.
A conventional method for correcting the bleaching has been done by multiplying the inverse of the ratio of intensity loss compared to a reference image frame. Estimation of bleach ratio could either be calculated directly 2-4 or by fitting exponential equation 5,6 . We call the first method simple ratio method and the second method exponential fitting method.
The simple ratio method seems to be the most widely used method. It has been applied for improving the deconvolution results in 3D stacks 7 and in the fluorescence recovery after the photobleaching (FRAP) technique. In FRAP literature, this method is called double normalization 2,3 . In ImageJ, this method was implemented as a macro by Jens Rietdorf and has been available since 2004 (https://www.embl.de/eamnet/html/bleach_correction. html).
In the exponential fitting method, it also compensates the loss by multiplying the inverse of the bleach ratio, but the method first fits the exponential equation to the bleaching curve and uses the fitted parameters for deriving the bleach ratio at each time point 8,9 . This method has been implemented as an ImageJ plugin and is available as a part of MBF ImageJ for Microscopy bundle (https://imagej.net/MBF_Plugin_Collection). A plugin PixBleach allows three different exponential decay models to be used for the correction (http://bigwww.epfl.ch/algorithms/ pixbleach/) 10 .
The histogram matching method is based on a strategy that is different from the above two methods and a novel algorithm introduced in this paper. Instead of correcting the fluorescence intensity based on the average intensity of the image at each time point, the histogram matching algorithm 11 is used to restore comparable intensity histogram distribution of each frame by taking the first frame as the reference.
We implemented an ImageJ plugin that allows bleach correction with a choice from these three different algorithms. In this article, we explain each of these algorithms in detail and compare their characteristics.

Implementation
In the case of the simple ratio and exponential fitting methods, the mean intensity of 3D stack from each time point was calculated to estimate the bleach ratio. For the histogram matching method, pixel intensity histogram was generated from the 3D stack at each time point.

Simple ratio method
We consider an i-th frame image I i = I i (x, y) in an image sequence and correct its loss of fluorescence emission. We assume that the mean intensity I is constant through the time-lapse sequence if not for the photobleaching. Then the ratio of the mean intensity of i-th frame i I to that of the first frame 0 I is the ratio of none-photobleached fluorophores in i-th frame. We could then estimate the true pixel intensity of the i-th frame by the equation below.  I I  I x y  I x y I  I I   −  = − − I b is the value of the background intensity. This value is estimated independently by measuring the none-fluorescence region within the image, or by measuring a blank image with all the image acquisition conditions being the same but without sample.

Exponential fitting method
The bleach correction tool included in MBF ImageJ bundle curve-fits total intensity value of each time frame with an exponential decay curve. This decay curve is then used to estimate the bleach ratio at each time point to calculate the true fluorescence intensity 12 . We implemented a similar capability for processing three-dimensional time series. In this case, the mean intensity of each time point, the average of 3D stack pixel values, was first fitted by an exponential equation to estimate the background intensity. again with the single exponential equation. Using the estimated parameters from this second fitting, which we now call them a′, b′, and c′, the ratio of bleaching was determined and then it's inverse was multiplied to the background-subtracted image.

Histogram matching method
Histogram matching algorithm modifies pixel values of an image to match its histogram shape to a reference image histogram 13 . We used the histogram of the first frame image H 0 (p) as a reference and matched the histogram of i-th image frame H i (p). p is the pixel value that is 0 ≤ p ≤ 255 in 8-bit image and is 0 ≤ p ≤ 65535 in 16-bit image. The cumulative distribution function of histogram CDF i (p) is used for the actual calculation.
Since we take the first frame of the time-lapse sequence as the reference, we use CDF 0 (p) as the reference CDF. We then match the rest of CDF to the reference CDF 0 (p) by

Operation
For running this plugin with Fiji, any entry-level laptop or desktop machine is sufficient. If one needs to work on a huge image stack, then one should make sure that the RAM has a capacity twice the size of the file size of the image stack. In case of Mac, Mac OS X 10.4 or higher is required to run Fiji.

Results
The sample image sequence was a time series of three-dimensional stacks 14,15 . The mean intensity showed an overall decrease by time, accompanied by repetitive small peaks (Figure 1, top-left). These small peaks corresponded to single time points, as each peak represented the spherical shape of the yeast cell. Mean intensity was low at the top slice, high at the cell equator, and then low again at the bottom slice.

Simple ratio method
To correct for the bleaching using the simple ratio method, we first determined background intensity. An arbitrary area outside the cell was selected and the mean intensity of the full sequence was measured. The mean intensity of the background was 68.3 ± 0.5. The sample sequence was then corrected for bleaching by the simple ratio method using this background intensity ( Figure 1, top-right). We further examined how the level of background intensity affects the correction results ( Figure 2). We used four different background values, 64, 68, 70, and 72. A slight difference in the background values caused a large difference in the resulting curves. When the background value was set large, the correction resulted in a curve with an increasing trend. When the value is set small the corrected curve showed a decreasing trend. When the measured background intensity 68 was used, drift was minimal.

Exponential fitting method
With the exponential fitting method, the estimation of background value is already a part of the algorithm. In the sample image, the estimated background intensity was 73. The image sequence was subtracted by this background value, fitted with a single exponential equation and then the bleaching was corrected (Figure 1, bottom-left). The mean intensity of the corrected image stack was mostly constant but with a slightly decreasing trend and also with small fluctuations. The decreasing trend was due to a fact that a single exponential equation was not perfectly fitting to the bleaching process seen in the sample. Since the exponential fitting assumes an idealized decrease in fluorescence by time, fluctuations present in the original image sequence were preserved in the corrected image sequence ( Figure 1, red curves).

Histogram matching method
The histogram matching method does not need background intensity estimation in its algorithm. For the purpose of comparing this method with the other two methods, the original stack was first subtracted by a constant background intensity value 68. The correction resulted in a stably constant mean intensity time series (Figure 1, bottom-right). It was the most stable result compared to the other two methods in terms of the fluctuation of intensity after the correction (Figure 1 red curves).

Discussion
All three methods do correct bleaching of fluorescence but the results of correction showed some difference. The simple ratio method has the capacity to correct time series with abrupt changes in intensity because the bleaching ratio is calculated for each time point. On the other hand, the quality of the correction is heavily dependent on the estimated value of background intensity. A small deviation of estimated value from the true background value causes wrong correction results (Figure 2).
The exponential fitting method has an assumption that the bleaching process follows a single exponential decay. It is known that in some cases bleaching time course is a double exponential decay 8 . In this respect, one must consider before using this method if there is any good reason to model the bleaching process of the sample as a single exponential decay. Otherwise, the goodness of fit needs to be evaluated for the proper use of this method. In addition, the exponential fitting method ignores small perturbations in the intensity such as abrupt changes in the emission of fluorescence. Such changes can be caused by small fluctuations in the power of the excitation light or in the slight variations in the timings of the shutter controlling the light path. The simple ratio method deals better with such changes. However, such a non-flexible nature of the exponential fitting method can become an advantage in some other occasions. For example, if the change in the intensity is due to the synthesis of GFP molecules by cell, the simple ratio method will wrongly correct such true increase and mask the biological event, but the exponential fitting method will achieve a better correction as long as the bleaching is known to be a single exponential decay.
The histogram matching method is robust when it is difficult to measure background intensity. This can happen when the whole image frame is filled with sample. For example, image data with packed cells in the image frame hinder the estimation of background intensity. Since histogram matching does not require the input of background intensity, the correction will be straight forward even with such cases. Moreover, this method is especially suited as a preprocessing for segmentation since strictly constant mean intensity in the corrected image sequence affords an optimal condition for segmenting objects.
The limitation of using the histogram matching method is that it assumes a stable distribution of fluorescence. If object under observation undergo changes in localization, e.g. signal changes from diffuse to spots during cell surface receptor internalization, or in shape, e.g. cell spreading, we expect that histogram shape will change as well. If we take an example of cell surface receptor internalization, the formation of spots by the aggregation of protein is expected to create a new peak in the high pixel values while decreasing the height of the existing peak in low pixel values of the histogram. Applying histogram matching to this time course would result in the wrong correction by forcing the histogram shape to become constant over time.

Conclusions
All thee methods correct bleaching, but have specific limitations of each. With the simple ratio method, background intensity should be accurately estimated. The exponential fitting method relies heavily on the model. With the histogram fitting method, object shape and pattern should be constant. For choosing an appropriate method, these limitations, and the known details of the observed biological event should be taken into account. In the future, the following features are planned to be added to the plugin: first, a fitting method with a double exponential equation; second, the background estimation method for the simple ratio method will be considered and will be added as a helper function; third, the current histogram matching method uses the first frame as the reference frame. Tolerance to changes in shape and pattern becomes higher if the reference frame is updated for every time frame so that the matching is done between neighboring time points.

Software availability
The Bleach correction plugin is contained within the download package of Fiji (https://imagej.net/Fiji) and can be used directly by launching the software package and accessing the menu item [Image > Adjust > Bleach Correction] without any further additional installation. For using it in ImageJ 1.× (ImageJ

BioImaging Faciliy & Division of Cell Biology, The Netherlands Cancer Institute, Amsterdam, The Netherlands
In this paper Kota Miura describes an ImageJ plugin for fluorescence photobleaching correction. This tool has been built into Fiji for a number of years already, and is no doubt used by many researchers. The plugin comprises three different methods to compensate for photobleaching in time-lapse movies (2D/3D). All methods alter the pixel intensities of the output image, in an attempt to stabilize the signal.
It is praiseworthy that the author takes the effort to describe such an established tool. The article thus also serves as a manual.
The article is well-written in a very straightforward way, without unnecessary details. The three methods are independently explained to a level that the reader can understand the workflow, but the descriptions don't go very deep. The discussion is to the point and mainly handles limitations of the implemented methods. The author provides adequate links to software, source code and sample data.
While reading the paper I had several questions and a few concerns, which I will discuss below. I realize that there are a lot of points; some are very technical and some can be seen as feature requests. I do not expect the author to address/implement them all.
All example data in the paper is taken from one sample image of a single cell. This image is acquired in 3D, but the bleach correction is performed, or at least shown, in a 'linearized' 2D stack. This 2D representation causes periodic 'peaks' in the time trace. Perhaps it is the author's intention to show that the tool works for 3D time-lapses, but these 'peaks' obscure the bleaching curve, are unnecessarily confusing the reader and hamper a thorough comparison of the methods. Luckily, the mean intensity of each time point is also provided (the red lines in figure 1 Because the output image is not in 32-bit floating point format, this operation also creates a positive bias in the background noise (negative numbers are set to 0). Even more so, the calculation itself seems to be done in 8-bit or 16-bit integer format. This partially underlies the following issue: ○ Figure 2 depicts corrected traces of a single cell surrounded by background, for varying input background estimation, using the Simple Ratio method. The figure seems to indicate that overestimating the background is much more severe than underestimating. However, the interesting part of the data is obviously the cell, and not the background. The corrected traces are mean intensities of the whole image. A flat trace can arise from the amplified signal of the bleaching cell alone (this is what we want), but also from an erroneously amplified background. With many more background pixels than 'cell pixels' this effect can be large even with moderate background increase. Exactly this turns out to be happening with the image stack presented here. I have reproduced the bleaching corrected data with the sample image used in the article, but now also separated the foreground (the cell) from background (using a simple Otsu threshold on a slice where the cell is in focus, halfway the time-lapse).  a, b, c... or a', b', c'?
Exponential method: The author correctly states that bleaching sometimes (more often than not?) follows a double exponential curve, and that "the goodness of fit needs to be evaluated for the proper use of this method." Unfortunately, the displayed R 2 in the plot is not the same as goodness of fit, and is really not valid for a nonlinear regression. (See for instance Spiess and Neumeyer (2010 1 ). A real goodness of fit, like Chi 2 , involves knowing the uncertainty of the signal. I am not asking for that, but for judging the goodness of fit by eye it would be instructive to also show a residuals plot. The article describes the Fiji plugin Bleach correction, which integrates three methods for correcting intensity loss by bleaching. One of the methods, histogram matching, is introduced and novel for addressing bleaching. All methods are well explained. Methods are compared to each other not only by performance on the given use case but also with other use scenarios in mind. I especially liked that specific biological examples were given (e.g. synthesis of GFP) and how the methods would compare in these scenarios.
It was no problem to test the plugin with the sample data linked in the text and the description.
Structure of the article. Methods and results could be fused, such that the explanation of a method goes directly along with presentation of the correction result on the sample data set. It should be referred to Fig 1 earlier in the text, e.g. Fig 1 top left to illustrate bleaching when introducing it in the very beginning. The chosen structure could be due to guidelines of the publisher and in this case it should be communicated to the publisher.
Other minor suggestions:

Introduction:
Include a sentence about ImageJ.

Methods: Simple ratio method:
We could then estimate the true pixel intensity I i c of ..
Include the information that the measurement is done per slice.

Histogram matching method:
Since the article has also educational character it could be good to include example cumulative histograms of a frame before and after correction using the histogram matching method.
We used the histogram of the first frame (per slice?) image H0 (p) as a reference..

Results:
The sample image sequence was a time series of three-dimensional stacks of a yeast cell.

Simple ratio method:
Add a brief definition of background to the introduction, and that, due to the different nature of source, it has a different bleaching behavior. I liked that the effect of the background intensity value was given much space in the article, including a figure.

Exponential fitting method:
In the sample image, the estimated background intensity was 73. : In the sample image, the background intensity estimated by fitting was 73 (to stress that estimation by fitting is different from the estimation by manually selecting a background region and measuring the intensity within the selection). You might want to show the fit together with the curve (as done in the plugin output). I have a few minor comments. The page numbers refer to the PDF version, as appear in the first version 1 of the manuscript.
Page 3, col.1 "before or after the experiment". Unclear, maybe add "imaging experiment", or simply "imaging" or "image acquisition" or "image data acquisition". The experiment, without context, is typically larger in scope and includes, for example, sample preparation.
○ "fluorescence loss is estimated and then the none[sic]-bleached condition is restored" is a bit too optimistic. Apart from the typo (non-bleached), I would argue something like "fluorescence loss is estimated and compensated for through image processing" would convey a more accurate assessment.
○ Page 3, col.2 "none-photobleached" => "non-photobleached" ○ "none-fluorescence" => "non-fluorescing"? or, better, "an empty region" ○ Page 4, Figure 1. The plot was confusing at first glance as it mixes two dimensions in one. I'd suggest using either a 3D surface (t, z vs Intensity) or an average/max projection in z (t, Intensity mean over z). The latter is probably sufficient even though it hides the Intensity variation in z.
Page 4, col. 1 There is no mention of the number of bins in the image. This is relevant for the following reason. CDF 0 can be inverted only if strictly monotonic (CDF 0 (a) < CDF 0 (b) i.i.f. a0 -1 will have two possible values c-1 and c. An easy fix would be to add a small constant to H everywhere (e.g. 1) after calculation.

Is the description of the software tool technically sound? Yes
Are sufficient details of the code, methods and analysis (if applicable) provided to allow replication of the software development and its use by others? Yes

Is sufficient information provided to allow interpretation of the expected output datasets and any results generated using the tool? Yes
Are the conclusions about the tool and its performance adequately supported by the findings presented in the article?
Yes alongside the graphs, for visual inspection.
The conclusion section clearly states the different fields of application/restriction use for each single method. A visual transcript, as a decision tree might help the end-user decide which method to try in her/his first attempts.
expertise to confirm that it is of an acceptable scientific standard.

Comments on this article
Version 1