Computer Science
Computer Science-2021.gz:143429
Computer Science
1276025bf4ed4b383cd1e1a3ba70575d52e5c743
Apoptosis quantification in tissue: development of a semi- automatic protocol and assessment of critical steps of image processing

Apoptosis is associated with numerous phenotypical characteristics, and is thus studied with many tools. In this study, we compared two broadly used apoptotic assays: TUNEL and staining with an antibody targeting the activated form of an effector caspase. To compare them, we developed a protocol based on commonly used tools such as filters, z-projection and thresholding. Even though it is commonly used in image-processing protocols, thresholding remains a recurring problem. Here we analyzed the impact of processing parameters and readout choice on the accuracy of apoptotic signal quantification. Our results show that TUNEL is quite robust, even if image processing parameters can allow or not to detect subtle differences of the apoptotic rate. On the contrary, images from anti-cleaved caspase staining are more sensitive to handle and proved to necessitate to be processed more carefully. We then developed an open source Fiji macro automatizing most steps of the image processing and quantification protocol. It is noteworthy that the field of application of this macro is wider than apoptosis as it can perfectly be used to treat and quantify other kind of images.


Introduction
Apoptosis is a programmed cell death characterized by caspases activation, subsequent degradation of cell components, including DNA fragmentation, and final phagocytosis of so called "apoptotic bodies" by surrounding cells or macrophages [1]. Importantly, apoptosis is not only critical for correct development of metazoan organisms, but also for their survival. Indeed, apoptosis failure is observed in many diseases including cancers. Therefore, it is widely studied and new actors are regularly identified. Apoptosis detection can be performed by multiple methods based on various features of apoptotic steps or regulators. Imaging of apoptosis in whole tissues can rely on a more limited number of methods. The first developed and best known of them is TUNEL (Terminal deoxynucleotidyl transferase dUTP Nick End Labeling) which is based on labeling of DNA 3' ends whose number increases during the DNA fragmentation step of apoptosis. However, TUNEL is costly, time consuming and also detects necrotic cells [2]. Alternatively, use of antibodies raised against cleaved ─ and thus activated ─ executioner caspases has proved to be more specific and convenient since immunodetection protocols are less time consuming as they include fewer steps than TUNEL. In mammals, the cleaved form of executioner caspase 3 is targeted [2]. In Drosophila melanogaster, the antibody used was raised against the executioner caspase Dcp-1 cleaved at Asp 216. This antibody was recently shown to actually detect the cleaved forms of both Dcp-1 and DrICE executioner caspases [3].
In this study, we wanted to compare TUNEL and cleaved caspase stainings the more objectively possible. To this end, we co-stained apoptotic wing imaginal discs (the larval tissue giving the adult wing) with TUNEL and anti-cleaved Dcp-1 antibody and addressed their sensitivity and requirements in terms of image processing.
As a single image can give a great diversity of information, the first step of image analysis consists in choosing a readout (i.e. which data is worth collecting). For example, protein quantity can be assessed by measuring staining intensity and tumor or bacterial colony growth can be followed by measuring stained area. In the case of apoptosis, the most commonly found readouts are the "number of apoptotic cells" or an "apoptotic index" that has various definitions according to the lab [4][5][6][7].
When the chosen readout is the "number of apoptotic cells", many studies use a manual counting, implying that an experimenter defines interesting spots and count them. Manual counting is reliable because the expert eyes of experimenters are able to distinguish the signal of interest from background noise better than any machine. This counting can be assisted by the "Cell Counter" plugin of ImageJ that records every experimenter's clicks. Even with this assistance, this approach remains time consuming and might involve estimation bias that can eventually raise ethical questions. This is why, whenever it is possible and whatever the readout, it is better to rely on automatized -or semi automatized-computer-based methods. However, machines do not have eyes trained to recognize specific signal from background. This discrimination is allowed by the image processing steps done prior to quantification in order to decrease background noise and amplify the signal of interest. This proper discrimination of foreground from background is called segmentation and defines the boundaries of the objects of interest. Thus, segmentation quality directly affects quantification accuracy.
Other software than ImageJ such as Imaris or Matlab display default functions for signal quantification [8,9] which usually comprise image processing to yield a rapid result. However, these programs are not open source and the methodology used to obtain the values is often hard to access, which means that users have only a limited control on their implementation. Moreover, they usually are computationally demanding and thus require powerful device to run. For all those reasons, many researchers prefer working on ImageJ / Fiji with which one can develop its own protocol for image processing and quantification.
Once the readout has been chosen, image quality has to be increased by getting rid of background noise and artefacts in order to improve segmentation. Images typically display three major kinds of defects: 1. general background noise; 2. isolated pixels with an aberrantly high intensity; 3. groups of pixels with aberrantly high intensities. Many functions are available to improve image quality on ImageJ but they often involve experimenter's appreciation. Unfortunately, the more the experimenter is involved, the harder it is to ensure that all the images have undergone the same process. However, most of the time, it is impossible to totally obliterate experimenter involvement.
We previously showed in the Drosophila model that overexpressing rbf1, the homolog of the human tumor suppressor RB1, induces apoptosis. This apoptosis requires the pro-apoptotic Bcl-2 family member Debcl, and involves caspases activation [10]. It can be visualized using TUNEL on rbf1 overexpressing wing imaginal discs [11]. In this study, we co-stained wing imaginal discs overexpressing rbf1 alone (vg > rbf1) or in the context of a debcl partial inactivation (vg > rbf1, debcl E26 ) with both TUNEL and anti-cleaved Dcp-1. We used these images to compare several methods of image processing and estimate their impact on the quantification of apoptosis for both assays. We then developed a semi-automatic protocol available as a free access Fiji macro called CASQITO (Computer Assisted Signal Quantification Including Threshold Options). This protocol enables, to process images of both labelings and quantify the number of apoptotic cells or the stained area. It is worth noting that our analysis and protocol can be relevant to quantify other types of staining outside the field of apoptosis.

Fly stocks
Flies were raised at 25°C on a standard medium. The UAS-rbf1 and vg-gal4 strains were generous gifts from Joel Silber (Institut Jacques Monod, Université de Paris, France). The debcl E26/E26 was obtained from the Bloomington Drosophila Stock Center (BL 27342) and we used a w 1118 fly stock as the reference strain.

Immunostaining and images acquisition
Third-instar larvae were dissected in 1X PBS pH 7.6 in order to remove every possible tissue except wing imaginal discs, then carcasses were fixed with 3.7 % formaldehyde in 1X PBS for 20 minutes at room temperature and washed three times for 10 min in PBST (1X PBS, 0.3 % Triton X-100). Discs still attached to cuticles were then saturated for 1h in PBST-BSA (1X PBS, 0.3 % Triton X-100, 2 % BSA) and dissected again to isolate wing imaginal discs which were then incubated overnight with 1:100 dilution of anti-cleaved Dcp-1 (Asp216, Cell Signaling Technology) at 4°C. The following day, after three washes in PBST, wing discs were incubated for two hours with anti-rabbit secondary antibody (1:400, Alexa-Fluor-612-conjugated goat anti-rabbit-IgG (H+L) antibody, Molecular Probes, Thermo Fisher Scientific) in PBST. Following three washes in PBST, TUNEL staining was performed according to manufacturer instructions (ApopTag Red In situ apoptosis detection kit, merck-millipore). Finally, wing discs were mounted in ProLong Diamond (Invitrogen) and images were acquired using a Leica SPE upright confocal microscope (Leica) at 568 nm for TUNEL and 612 nm for anti-cleaved Dcp-1 stainings. Image analysis was done exclusively on Fiji, the exact same zone selection was studied for both assays 3. Results

Choosing a readout according to the biological question
The readout is the data used to translate the intensity of the biological effect in numbers. Therefore, the chosen readout should be coherent with the biological question and the tool used to study it. For instance, intensity can be measured to assess the amount of a stained component. In the study of apoptosis rate, whatever the assay used, a cell is apoptotic or not. Thus, quantifying the intensity of the staining, even if it can somehow make sense, does not seem the best option for accurate quantification of apoptosis. Conversely, as long as apoptotic cells can be separated from each other (low apoptosis rate, widespread pattern or intracellular discrete staining), counting the number of objects equals counting the number of apoptotic cells, which constitutes a valid readout. In case this readout cannot be used, another valuable readout is the stained area. This can be used as a readout per se or can be used as a primary data and further treated to get an apoptotic index or score. If so, the stained area can be divided per the number of cells (obtained by plasma membrane or nuclear co staining), the area of interest (surface of a cellular clone or of the tissue section). These numbers do not indicate the actual number of apoptotic cells but this number is rarely necessary and those readouts satisfy the need to have a quantification precise enough to compare different samples.
Here, we used TUNEL and anti-cleaved Dcp 1 to detect apoptosis. These highlight different features of apoptosis as TUNEL labels fragmented DNA in the nuclei while anti-cleaved Dcp 1 staining is cytosolic. As rbf1 overexpression is a potent apoptosis inducer in the wing imaginal discs, the probability to have clusters of adjacent apoptotic cells is rather high. This can eventually become problematic for accurate quantification. Indeed, when adjacent cells are apoptotic, TUNEL labeling is expected to remain punctiform as nuclei remain spaced by cytoplasms (Figure 1 (a) and (c)). On the contrary, with anti-cleaved Dcp 1 staining, it is expected that such adjacent apoptotic cells become indistinguishable from each other and thus appear as a single object (Figure 1   Therefore, for the count of apoptotic cells, these clusters of labeled cells are not expected to alter the quantification for TUNEL while they may cause an underestimation of the number of apoptotic cells with anti-cleaved Dcp-1 staining. The extent of this underestimation is difficult to anticipate as it depends on many parameters. Still, this underestimation surely increases with the apoptotic rate -as the probability to have clusters of apoptotic cells increases -which could lead to an artificial flattening of the difference of apoptosis rate that may exist between two conditions. As for the area readout, the size of the wing imaginal disc cells (and their nucleus) being homogenous, the stained area indirectly reflects the number of apoptotic cells without being impacted by their relative localization. In the end, counting cells seems, at least at first sight, a more precise, because more direct, readout of apoptosis than area. However, this readout might be altered by apoptotic cells clusters. As it is not possible to anticipate how these clusters will affect the quantification in our experimental set up, we chose to use both count and area readouts. On Fiji, these two readouts can be obtained using the "Analyze Particles" function, which only works on binary 2D images. This means that our image processing protocol should include both z projection and binarization using a threshold, these two treatments being compatible with our set up. Indeed, our tissue of interest is a monolayer, then z projection should not affect quantification. Besides, because a cell is apoptotic or not, our readouts do not depend on signal intensity and binarization by itself should not affect the quantification.

Designing an image processing protocol
Steps of image processing directly depend on the chosen readout. In order to get both the number of objects and the stained area, our image processing protocol is based on three major steps: 1. Background noise reduction, 2. Compression of our 3D images into 2D by a z-axis projection, 3. Thresholding. Those steps allow appropriate segmentation required for relevant quantification by the "Analyze Particles". Importantly, on ImageJ, there are many ways to minimize background signal, 6 ways of compressing a 3D image in 2D and 17 ways of determining a threshold, resulting in countless combinations of possible image processing. In this study, we investigated the weight of these parameters on signal segmentation to end up with an optimized and unbiased protocol for apoptosis quantification.

2. 1. Median Filter and size limitation efficiently reduce artefacts
When quantification is automatized, definition of the signal of interest by segmentation is even more critical. Indeed, bad segmentation can lead to quantification of unreal objects and thus gives useless results. To do so, signal of interest boundaries have to be better defined while background noise have to be decreased. General background noise can be minimized in many ways depending on the kind of images, the readout wanted and the defaults faced. In our case, mandatory use of a threshold would blacken every low intensity pixel responsible of general background noise. However, if binarization of the image efficiently removes diffuse low intensity background noise, it is not sufficient to erase artifactual pixels with aberrantly high intensity, i.e. which intensity is higher than the threshold value. Fortunately, isolated aberrantly high pixels can be dealt with filters. Filters are matrix operations that re-calculate a pixel intensity value based on itself and its neighbors. The two mainly used filters are "Mean Filter" and "Median Filter" (Figure 2A). . In order to make intensity differences more visible, each boxes background color corresponds to the double of each pixel intensity value in greyscale. B: (a) TUNEL-labeled wing imaginal disc image after a Max Intensity z-projection. (b) magnification of (a) after binarization using a manually determined threshold. (c) same as (b) but the image was applied a "Median Filter" with a radius of 1 before Max Intensity z-projection. White bars correspond to 10 µm. Arrows with circled numbers 1, 2 and 3 target areas of interest.
A "Mean Filter" with a radius of 1 gives a pixel an intensity value that corresponds to the mean of its value and those of its direct neighbors (Figure 2A (b)). Hence, the value of an isolated pixel with aberrantly high intensity is attenuated by the intensity value of its neighbors. However, as extreme values impact mean calculation, every neighboring pixel is affected by the isolated aberrantly high pixel and their intensity is artificially increased (Figure 2A (b)). By contrast, a "Median Filter" with a radius of 1 gives a pixel an intensity value corresponding to the median of its value and those of its direct neighbors (Figure 2A (c)) which is expected to be much closer to local intensity value. Besides, as extreme values effect on median calculation is low, high intensity isolated pixel impact on its neighbors is negligible. All in all, the "Mean Filter" tends to spread an aberrantly high value whereas a "Median Filter" tends to confine it. This effect is illustrated in Figure 2B, where the "1" arrows of the "No filter" panel shows typical isolated aberrantly high pixels that are efficiently erased by a "Median Filter" with a radius of 1 ( Figure 2B compare (b) and (c)). Aside from this benefit, the area pointed by the "2" arrow exemplifies the median filter ability to preserve edges of an object. Indeed, on the illustration, human eyes easily detect that the "2" arrow targets a marked cell ( Figure 2B (a)). However, this object is heterogeneous: in a restricted space, it contains few pixels of high intensity and many pixels of low intensity (i.e. below the chosen threshold). Without any filter, only high intensity pixels are kept after the thresholding step thereby fragmenting this object in several small groups of pixels ( Figure 2B (b)). Thus, with no further treatment, multiple objects will be counted in this area, which does not reflect reality. However, as these high intensity pixels are close to each other, the "Median Filter" with a radius of 1 homogenizes intensity values within this object. This allows its reconstruction and gives a segmentation consistent with reality ( Figure 2B (c)). Thus, when the "Median Filter" with a radius of 1 is applied, the "number of objects" decreases only to be closer to what a human eye would count.
Although the "Median Filter" with a radius of 1 efficiently reduces the number of artifactual objects by erasing isolated high pixels, the issue of groups of pixels with an intensity higher than the threshold value remains. A way to eliminate most of those artifacts is to limit our analysis to objects with a size consistent with the smallest biological object of interest. In our case, this smallest biological object is TUNEL-labeled nucleus, we assessed their size on a few random images and thus set a size limit at 2 µm. Importantly, this "> 2 µm size limitation" fits our data but should not be taken as a default value and must be adapted for other kinds of signals or cell types. As the "Analyze Particles" function records the size of every object, this filtering can be done after quantification. This function directly proposes to define a size range of object to quantify which we used. In Figure 2B, the "> 2 µm limitation" eliminates artifactual object pointed by the "3" arrow as well as individual pixels such as those pointed by the "1" arrows. It thus appears very powerful to "clean" the image. However, as efficient as the size limitation may be, it cannot replace "Median filter". Indeed, as already explained, in the absence of a "Median Filter", the cell indicated by the "2" arrow in Figure 2B gets fragmented in several small groups of pixels, each one being smaller than 2 µm ( Figure 2B (b)). Thus, without the "Median Filter", these pixels are eliminated by the "> 2 µm limitation" and the actually labeled cell indicated by the "2" arrow is not included in the quantification of the apoptotic signal. Here, reconstruction of the object by the "Median Filter" prevents its elimination by the "> 2 µm limitation" (Figure 2B (c)). In the end, combination of a "Median Filter" with a radius of 1 and "> 2 µm limitation" allows a better segmentation and a more accurate quantification.

Max Intensity z-projection improves contrast
Confocal microscopy gives the possibility to capture objects in 3D. However, image processing often requires transforming volumes into 2D images by compressing the z-axis. In our case, the "Analyze Particles" function used to quantify the signal of interest The z-axis presented in "A" shows only one illuminated pixel on slice #2 and this pixel is included in the object. The "B" z-axis shows numerous pixels highly illuminated, comprised in the object. The "3" z-axis presents an artifactually illuminated pixel on slice #2 that is not comprised in the object. Comparison of the projection methods shows that AI projection decreases the importance of the artifactual pixel of the "C" z-axis while the MI projection, increases its weight. However, in the processing protocol, the preceding use of a "Median Filter" with a radius of 1 eliminates most of those artefacts that are thus not present anymore at the projection step. Conversely, the AI projection of the "A" z-axis leads to loss of signal even if it is part of the object. Furthermore, with the AI projection, the contrast between object and background is little so the range for the appropriate threshold value is limited (Figure 3). After a MI projection, contrast is enhanced and thus threshold determination is easier for the experimenter, which helps limiting the experimenter bias. This is particularly important for signals with low contrast such as TUNEL. In our case, these two projection methods do not end up in drastically different results but, all in all, MI projection presents more benefits than AI projection.

Use of custom manual thresholds gives the best segmentation for relevant quantification
The "Analyze Particles" function used for quantification requires the image to be binary. The transition from a greyscale image to a black and white image involves the setting of a threshold that defines an intensity value above which a pixel is turned to white and under which a pixel is turned to black. Ideally, this value should enable to get an image where white only corresponds to the signal of interest. Thresholding is the last step of segmentation and finally defines the objects of interest, which is critical for accurate quantification. Therefore, among the steps of image processing, it is the one that has the most dramatic effect on quantification accuracy so we dedicated a particular attention to the threshold determination method.
Threshold can be automatically set by algorithms that analyze specific features of the image intensity histogram to determine a threshold value using either simple indicators such as the mean, maximal or minimal intensity values, or more complex formulas. Hence, algorithms appear as an unbiased method to obtain a specific threshold value per image. We thus wondered whether any of the 16 thresholding algorithms available on Im-ageJ / Fiji could be used to determine a threshold capable to properly segment apoptotic signal on our images. Using some randomly chosen images, we visually checked if these algorithms could provide a threshold value allowing a relevant segmentation, i.e. consistent with apoptotic signal. Most algorithms did not pass the visual inspection step as they yielded unrealistic segmentation either by ignoring a great portion of the signal or by including artifactual signal. However, two of them, Otsu and Moments, seemed capable to discriminate actual apoptotic staining from background. We then performed a more detailed analysis of the threshold values obtained with these algorithms by comparing them to the ones obtained by experimenters. To this end, for the 28 images of the vg > rbf1 genotype, three experimenters determined the threshold to use for each staining (TUNEL or anti-cleaved Dcp-1) by eye and in triplicate (see Supplementary Figure S1). Thresholds for these images were also determined using the 16 algorithms. As expected, algorithms inducing obvious unrealistic segmentation of the apoptotic signal yielded threshold values very far from the range of the ones determined by experimenters (Figure 4   From this global analysis, it could seem that algorithms can be as good as human eye for threshold determination (compare for instance Otsu and Exp in Figure 4 (b)). However, visually, Otsu capability to determine a relevant threshold seemed irregular. We thus further digged into these data and compared the threshold values obtained not globally but for each image. As shown in Figure 4 (c) and (d), the values obtained with Otsu, if they tend to be roughly the same on average than the ones obtained by experimenters, are actually most of the time out of the range of experimenters values. This is particularly striking for images obtained from TUNEL (Figure 4 (c)) as Otsu's determined values are far higher or lower that the ones obtained by any experimenter. This would not be an issue as long as the threshold values obtained still allow a realistic segmentation of apoptotic signal and subsequent relevant quantification. However, such deviation of the threshold value results in an inappropriate segmentation (compare Figure 4 (e) and (f)), that necessarily ends in a biased or most likely totally wrong quantification. When it comes to Moments, it provides threshold values that are usually lower than experimenters' ones ( Figure 4 (a) and (b)), which means that using this algorithm tends to include some background noise to the quantification. The question resides then in determining whether this amount of background noise is important enough to alter quantification. In the case of TUNEL staining (Figure 4 (c)), values are quite low so it certainly affects quantification rather importantly. By contrast, when images come from anti-cleaved Dcp-1 staining, the threshold values given by Moments are much closer to the ones obtained by experimenters, they actually seem very similar to the lowest values determined by experimenters (Figure 4 (d)). Therefore, one could assume that the variability of threshold values between Moments and an experimenter is comparable to the one that exists between experimenters. We tested this by comparing the relationship between the two most distant experimenters' batch of measurements to the one between Moments and its closest dataset. As shown in Supplementary Figure S2B, if two experimenters will not determine exactly the same value for the threshold, their evaluations remain consistent with each other (p = 10 -5 and R 2 = 0.53 for the most distant measurements), the difference can be more or less described as a given experimenter tending to set thresholds always lower than the other. This is a systematic error that should affect quantification only moderately. On the contrary, threshold values determined by Moments are not consistent with the values of experimenters (p = 0.8 for the closest in Supplementary Figure S2B). This indicates that Moments can set a low threshold value when experimenters would have all chosen a higher one but it is not always the case, and most of all, the extent of this underestimation (i.e. the range of the difference between experimenters and Moments threshold values) is variable. This is more problematic as the extent of background incorporation in the quantification will then vary and might alter quantification relevance.
Contrary to algorithms, manual determination of the threshold values appears quite robust. Indeed, comparison of manually determined threshold values for individual images shows a low variability and a good reproducibility both between several determinations of a given experimenter and between experimenters (Supplementary Figure S1). As all images from an experiment are acquired identically, with the same microscope settings, originate from samples treated with the same solutions, at the same time, theoretically, the appropriate threshold value is expected to be the same for all of them. Moreover, using the same unique threshold value for all images can be considered as more objective and unbiased.
However, manually determined threshold values display some variability (Figure 4, Exp). This can justify using a distinct individual threshold for each image as it might enable a more accurate segmentation and subsequent quantification. In order to assess to which extent these two thresholding methods can affect quantification and detection of our biological effect, we tried both ( Figure 5A).
In the Manual condition, each image was binarized using its own manually determined threshold value (determined by experimenter 1, measure 3). From these individual threshold values, we calculated the median value per genotype and then the median of these medians. This last value was used as the unique threshold value to binarize all images in the Total condition (1128 for TUNEL and 777 for anti-cleaved Dcp-1). We decided to use the median of the medians per genotype rather than the global median (calculated from the whole of the images independently of their genotype) to avoid giving more weight to a genotype (that might have a larger headcount for instance).  a) and (b)) or anti-cleaved Dcp-1 ((c) and (d)) for vg > rbf1 ((a) and (c)) and vg > rbf1, debcl E26 ((b) and (d)) genotypes. White bar corresponds to a 50 µm scale. "Count" panel presents quantification of TUNEL (e) and anti-cleaved Dcp-1 (f) signal based on the counting of the number of objects according to the different thresholding methods (Manual and Total) for vg > rbf1 (in pink) and vg > rbf1, debcl E26 (in blue). "Area" panel presents quantification of TUNEL (g) and anti-cleaved Dcp-1 (h) signal based on the number of white pixels (stained area) according to the different thresholding methods (Manual, or Total) for vg > rbf1 (in pink) and vg > rbf1, debcl E26 (in blue). pvalues displayed above results were obtained using Wilcoxon tests.
For TUNEL stained images, both thresholding approaches give a workable quantification both for count and area readouts ( Figure 5B (e) and (g)). Indeed, a significant decrease of apoptosis between the two genotypes is detected in all cases. However, even if Total or Manual thresholding method enable to detect the biological effect, we noticed that quantification is still somehow altered when a unique threshold is used (Supplementary Figure S3).
As for anti-cleaved Dcp-1 stained images, the first observation we can do is that the count readout was not usable. Indeed, we knew that apoptotic cells clusters might alter quantification as these clusters might be considered as a single object. Moreover, such underestimation is enhanced when the apoptosis rate rises, eventually leading to the flattening of the difference between two samples. However, we chose to keep this readout, as it was not possible to anticipate the extent of this phenomenon in our samples. As shown in Figure 5 (f), with this readout, the difference in the apoptosis rate between the two genotypes becomes undetectable. This indicates that the level of apoptosis induced by rbf1 generates apoptotic cells clusters frequently enough to significantly alter quantification, and this, whatever the thresholding method, thus prohibiting the use of the count readout. By contrast, when anti-cleaved Dcp-1 staining is quantified using the area readout, the difference between the two thresholding methods (unique versus individual thresholds) becomes obvious. As shown Figure 5 (h), when a unique threshold value is used for all images (Total), the difference of apoptosis rate between the two genotypes is barely detectable (p = 0.033). Moreover, extreme values compatible with an overestimation due to inadequate segmentation are seen (see highest values for vg > rbf1, debcl E26 genotype in Figure 5 (h) and also Supplementary Figure S3). On the contrary, the use of individual specific threshold values (Manual) enables to readily detect the difference of apoptosis rate between the two genotypes (p = 1.5 10 -5 ).
In the end, this analysis shows that, in our case, using an individual threshold per image is more adequate and turns out to be the safest option for accurate segmentation and thus, relevant quantification.

Quantifications of TUNEL or anti-cleaved Dcp-1 stainings do not have the same requirements
TUNEL and anti-cleaved caspase stainings are widely used to assess the level of apoptosis in tissues. However, depending on the experimental set-up, the quantification step may become tricky.
TUNEL appears as a quite robust apoptosis detection technique. Indeed, it allows to quantify apoptosis and to detect differences in apoptosis rate whatever the thresholding approach, and with both count and area readouts ( Figure 5 (e) and (g)). This was not totally expected since in our images, there was not a strong contrast between the apoptotic signal and the background (Figure 5 (a) and (b)). However, as previously mentioned, TUNEL assay is costly, time consuming and lacks specificity as it also detects necrotic cells.
On the contrary, using antibodies against cleaved caspase(s) is considered as a more specific and convenient staining of apoptotic cells. By contrast with TUNEL which labels nuclei, caspase staining covers the whole volume of the cell, raising the issue of adjacent apoptotic cells when the readout is counting cells. Indeed, counting the number of apoptotic cells stained with anti-cleaved caspase antibody is perfectly possible as long as apoptotic cells are sufficiently separated from each other. In our experimental set-up it appeared that the apoptosis rate was too high to prevent the underestimation of the signal due to fusion in a single object of clustered apoptotic cells. When the stained area was measured, it revealed that images from anti-cleaved Dcp-1 stainings should be carefully processed because determination of the threshold value to use for binarization is particularly important. Indeed, even if the difference of apoptotic rate between the two genotypes was known and easily seen by eye ( Figure 5, compare (c) and (d)), its detection after quantification was not obvious. Actually, the decrease in apoptosis between the two genotypes is barely detected when a unique threshold value is used for the segmentation of all images ( Figure 5 (h)) whereas using a specific threshold value for each image enables to see it. Therefore, anti-cleaved Dcp-1 staining quantification is more affected by image processing than TUNEL and should be handled more carefully. Figure 6. Fiji macro workflow. Giving our image processing protocol, the macro is divided in two major parts. Part 1 (left panel) is dedicated to manual threshold determination while Part 2 (right panel) is dedicated to zone selection and quantification after application of the previously determined threshold. In Part 1, images are opened, treated according to parameters set in the window presented Figure 7 and presented to user for threshold determination. Once every image of a .lif file has been treated in Part 1, threshold results are recapitulated before starting Part 2. In Part 2, images are opened again and treated as in Part 1, chosen threshold is applied and the resulting image is presented to user for zone selection before quantification. When every relevant images have been treated, quantification results are superficially analyzed to yield mean, standard deviation, min and max.

Macro explanation
Once a .lif file has been chosen, the window presented Figure 7 opens to set a few parameters. Figure 7: Macro Part 1. The macro starts by collecting few parameters. Importantly, the chosen "Condition name" will end up in the name of the files associated with the .lif file. The default background noise reduction method is a "Median Filter" with a radius of 1 but user can apply its own protocol before and / or after Z projection if needed. As a thresholded image is in black and white, it can hardly be used to define region of interest when a zone selection is needed. Therefore, the macro offers the possibility to define region of interest (i.e. vestigial domain in our case) on another channel or an unthresholded version of channel of interest. If the whole image has to be analyzed, a "No selection needed" option is available. Concerning threshold determination, "A. Using an algorithm" and "B. Manual determination" options will lead to application of an individual threshold value per image. If use of a unique, representative, threshold value to treat every image of a .lif file is wanted, user can make its choice after treating every image (option "C"), 50 % of the images (option "D"), 25 % of the images (option "E") or 10 % of the images (option "F"). If manual determination has already been done, user can skip Part 1 with option "G". Once chosen, those parameters are recorded in the Log window which content is ultimately saved in a .txt file.
These parameters aim at setting the method used for background noise reduction and z-projection as well as determining channels of interest. It should be noted that channels are automatically assigned a number that increments in ascending order of the wavelength starting from the shortest wavelength used. Once these parameters are set, the macro opens the first image and checks if it is an actual stack and not a 2D image. In the latter case, the 2D image is skipped without being analyzed and the macro records it has ignored this image in the Log window ("Caution series X is not a stack") and directly opens next image. When the opened image is a stack, it is processed according to the parameters previously chosen for background noise reduction and z-projection. In our case, this means a "Median Filter" with a radius of 1 is applied and z-projection is done by the Max Intensity method. Result of this z-projection opens twice, giving one image to work on and a second to check in real time if the changes induced are consistent with reality. Just before allowing the user to determine a threshold value, a checkpoint asks for confirmation to treat this image. Indeed, acquisition might have been done on another channel than the one of interest thus, huge artifacts on channel of interest might reveal themselves only at this point. Discarding such image from the analysis is done by clicking the "Do not take series X into account" option. The macro records this image has been ignored in the Log window ("Series X has been excluded from threshold determination") and directly opens next image. At this step, four windows are open to allow threshold determination: the two images obtained after z-projection to monitor the effect of the threshold and two system windows: a "Threshold" window that contains the cursor used to set the threshold value and a "Threshold Selection" window whose "OK" button clicking is required once the threshold has been determined. Importantly, if the "OK" button of the "Threshold Selection" window must be clicked for the macro to continue, buttons of the "Threshold" window should not be clicked. As soon as the "OK" button of the "Threshold Selection" window is clicked, the chosen threshold value is recorded in the Log window and the .xls table and next series is opened. Once all images have been processed, Part 1 is over and a histogram showing the threshold values distribution is displayed. When ready, user have to click the "OK" button of the "End of Part 1" window to start Part 2. Importantly, Part 2 does not necessarily have to be done right after Part 1. As presented in the Part 1 settings choice window (Figure 7), a possible option of the last question is "G. Part 1 already done". This option will lead directly to Part 2. In this case, user can either define a unique threshold value to treat every image of a .lif file or treat each image with its own threshold value after loading the corresponding Excel table containing threshold values for this .lif file obtained in Part 1. This option allows users wanting to apply the same threshold value for all of their conditions to assess the best threshold value for all their conditions in a first time and perform quantification in a second time.
Part 2 begins with the opening of the window presented Figure 8, which offers the possibility to use the threshold value previously determined in Part 1 and to apply additional limitations to quantification such as a size limitation. The one presented here corresponds to option "B. Manual determination". If option "A. Using an algorithm" was chosen, the window would start by asking which algorithm should be used. If any of options "C", "D", "E" or "F" were chosen, the window would start by asking which unique threshold value should be used. The following parameters (how to retrieve the results, objects' size, objects' circularity and range of series to analyze) are always present. Once chosen, those parameters are recorded in the Log window which content is ultimately saved in a .txt file.
Furthermore, depending on the study, the data needed can be more or less detailed. Here, a global count of the number of objects or total stained area per image was sufficient to be conclusive. Therefore, only a summarized compilation of data was needed. When checked, this option generates an .xls file where each line corresponds to a treated image and its summarized result (i.e. total number of objects or total stained area in the ROI). However, one might need to compare objects size within an image and thus need more detailed data where specifications of every object are recorded. When checked, this option generates an .xls file where each line corresponds to an object (objects of all images are in the same table). Both options are available on the macro and user can choose one or both of them (in which case two .xls files will be generated). Once these settings are done, images are processed as in Part 1 except that, this time, the chosen threshold is applied. After thresholding, Part 2 offers the possibility to limit quantification to a region of interest. If so, the channel needed to draw the selection is opened. Then, when the "OK" button of the "Zone Selection" window is clicked, the macro goes on with quantification and results are recorded in an .xls file. If the settings chosen in Part 1 end up associating a threshold determined by an algorithm with a zone selection, the zone selection is done before running the algorithm. Indeed, to determine a threshold value, algorithms take the whole image into account, which can become a problem if there is a highly illuminated artifact outside of the region of interest. To allow a wide range of application for this macro, we have chosen to ask the "Analyze Particles" function to quantify all the possible readouts.

Conclusions
Apoptosis quantification in a tissue is usually indirect as it generally relies on imaging techniques. There are many ways to analyze an image and once a readout has been chosen, many processing protocols are possible. Here, we describe a semi-automatic protocol running on Fiji for quantification of apoptosis on Drosophila wing imaginal discs after TUNEL or activated-caspase labelings. During the development of this protocol, we paid particular attention to the weight of specific steps to obtain a realistic segmentation, which underlies an accurate quantification. As in many image processing protocols, determination of the threshold for binarization turned out to be a critical step. In our case, none of the algorithms available in Fiji was satisfying to determine relevant thresholds.
We also considered using the same threshold value to treat several images but, in the end,, the best option for our data, was to use a specific threshold manually determined for each image. Indeed, this method proved to carry out a proper segmentation for all images resulting in valid quantification and subsequent detection of biological effects. Even if one could be concerned about the bias that might be induced by this approach, the bias is in fact limited as we observed that threshold values obtained by experimenters are actually very consistent both for a given experimenter and between experimenters. Moreover, an appropriate processing of the images can facilitate this determination of a threshold value. In this sense, association of a "Median Filter" with a radius of 1 and a Max Intensity z-projection proved to be very efficient. It would also be interesting to try the "Sum Slices" projection that adds up all pixels intensity of a z-axis which should enhance contrast even more (particularly after a median filter) and thus facilitate threshold determination.
The protocol presented here should not affect other readouts available in the "Analyze Particles" function such as: bounding rectangle, shape descriptors, centroid, perimeter, Feret's diameter or stack position. Furthermore, we kept the options implemented in CASQITO macro rather open offering a possible use of this tool for a great variety of readouts, stainings and biological questions.