Root Mean Square Minimum Distance as a Quality Metric for Stochastic Optical Localization Nanoscopy Images

A localization algorithm in stochastic optical localization nanoscopy plays an important role in obtaining a high-quality image. A universal and objective metric is crucial and necessary to evaluate qualities of nanoscopy images and performances of localization algorithms. In this paper, we propose root mean square minimum distance (RMSMD) as a quality metric for localization nanoscopy images. RMSMD measures an average, local, and mutual fitness between two sets of points. Its properties common to a distance metric as well as unique to itself are presented. The ambiguity, discontinuity, and inappropriateness of the metrics of accuracy, precision, recall, and Jaccard index, which are currently used in the literature, are analyzed. A numerical example demonstrates the advantages of RMSMD over the four existing metrics that fail to distinguish qualities of different nanoscopy images in certain conditions. The unbiased Gaussian estimator that achieves the Fisher information and Cramer-Rao lower bound (CRLB) of a single data frame is proposed to benchmark the quality of localization nanoscopy images and the performance of localization algorithms. The information-achieving estimator is simulated in an example and the result demonstrates the superior sensitivity of RMSMD over the other four metrics. As a universal and objective metric, RMSMD can be broadly employed in various applications to measure the mutual fitness of two sets of points.

A localization algorithm or an estimator that estimates emitter locations plays an important role in obtaining a high-quality image in stochastic optical localization nanoscopy. A number of localization algorithms and software packages have been developed in the literature. It is imperative to identify the good algorithms, which yield a high quality of nanoscopy images, and to advance research and development of localization algorithms. To achieve this goal, by challenging 2D synthetic nanoscopy data, performances of thirty two localization software packages were recently evaluated 1 . The challenge has been advanced to focus on 3D imaging 2 and become an open public online challenge that up to now has drawn eighty four participant packages 3 . A quality metric for localization nanoscopy images is crucial and necessary in research and development of localization algorithms. However, a localization nanoscopy image is substantially different from a conventional image. The former consists of a set of points defined over an n-dimensional real space while the latter consists of an array of pixels with a finite number of intensities. A quality metric for the conventional images usually fails to evaluate the quality of the localization nanoscopy images. To emulate a quality metric for the conventional images, Fourier ring correlation (FRC) has been proposed to evaluate resolution of localization nanoscopy images 4,5 . Cramer Rao lower bound (CRLB) 6,7 can be employed to evaluate the minimum variance between an estimated emitter location and its true location for all unbiased localization algorithms. Nevertheless, these metrics are unable to manifest the quality of a localization nanoscopy image as an estimate of a set of emitter locations. A universal and objective metric, which is crucial and necessary to evaluate qualities of nanoscopy images and performances of localization algorithms, has not yet been established in the field 2 . In the current research, accuracy, precision, recall, and Jaccard index (JAC) are used as quality metrics for nanoscopy images 1,3,8,9 . These metrics depend on the full-width half-maximum (FWHM) of the point spread function (PSF) in an optical system and therefore are not universal and somewhat subjective, and in certain conditions fail to distinguish qualities of different nanoscopy images. In this paper, we propose root mean square minimum distance (RMSMD) as a quality metric for a localization nanoscopy image when the set of emitter locations is known. RMSMD is defined as root mean square of minimum distances from a point of one set to another set and vice versus. RMSMD measures the average, local, and mutual fitness between two sets of points. The smaller the RMSMD, the better the fitness or the higher quality of one set in description of another. The Voronoi cell of a point can be employed to express RMSMD. RMSMD has several properties that are common to a distance metric as RMSMD is. The kernel points of the two sets play an important role in RMSMD since removal of non-kernel points can usually reduce RMSMD. A numerical example demonstrates that RMSMD can properly measure the quality of a sequence of localization nanoscopy images while the metrics of accuracy, precision, recall, and JAC fail to distinguish the qualities of different nanoscopy images in certain conditions. The ambiguity, discontinuity, and inappropriateness of the metrics of accuracy, precision, recall, and JAC are analyzed. A simulation of experiment that an unbiased Gaussian estimator achieving the Fisher information (therefore achieving CRLB) localizes emitter locations frame by frame independently is carried out. The simulation result demonstrates the superior sensitivity of RMSMD to a quality change over the other four metrics. The RMSMD and visual quality of the information-achieving localization nanoscopy image can be used as a benchmark for all localization nanoscopy images and algorithms. As a universal and objective metric, RMSMD can be broadly employed in applications to measure the fitness between two sets of points.

Method
Definitions. Consider two sets of n-dimensional points S = {s i , i = 1, …, M} and X = {x i , i = 1, …, N} where s i 's or x i 's are unnecessary to be distinct. The mean square minimum distance (MSMD) between X and S is defined as where |·| is the number of elements in a set and ⋅   is the l 2 norm or the Euclidean distance between two points. Then their root mean square minimum distance (RMSMD) is D(X, S). If a point in one set has the same minimum distance to multiple points in the other set, one of the multiple points is arbitrarily chosen and the RMSMD does not depend on the choice.
In the definition, − ∈ x s min x X 2 for s ∈ S is the minimum square distance from the point s to the set X, and − ∈ s x min s S 2 for x ∈ X is the minimum square distance from the point x to the set S. Therefore, D 2 (X, S) is the average minimum square distance from a point in X to S and from a point in S to X. Because of this, D(X, S) evaluates how well the two sets X and S averagely, locally, and mutually fit to each other. The smaller the RMSMD, the better the fitness between the two sets. Since D(X, S) depends only on X and S regardless of how they are produced, RMSMD is a universal and objective metric and can be extensively applied in various applications.
It is worth noting that both summation terms in the definition of RMSMD are necessary; otherwise it may result in an awkward measurement of fitness between X and S. For example, suppose that only the second summation term in Eq. (1) were considered in the definition of RMSMD. When M ≥ 2, N = 1, and x 1 = s 1 ≠ s 2 the second term is equal to zero, which would imply a perfect fitness. However, X and S could be actually significantly different in the case that s 2 is far away from s 1 . Such a difference is taken into account by the first summation term in Eq. (1).
The Voronoi cell of s i ∈ S is defined by that contains all the points in  n , whose distance to s i is not greater than to any other s j for j ≠ i. Any two adjacent Voronoi cells V(s i ) and V(s j ) are separated by the hyperplane that contains the middle point of and is orthogonal to the line segment connecting s i and s j . The n-dimensional space  n is partitioned by the M Voronoi cells V(s i ). Similarly, the Voronoi cell V(x i ) for x i ∈ X is defined. The n-dimensional space  n is also partitioned by the N Voronoi cells V(x i ). In terms of the Voronoi cells, the MSMD can be written as We define a kernel point and a kernel set as follows. x * is the point in X that has the minimum distance to a point s * in S, i.e., = − ∈ ⁎ ⁎ x x s arg min x X ; meanwhile, the point s * is also the point in S that has the minimum distance to x * , i.e., = − ∈ ⁎ ⁎ s s x arg min s S . Then x * is said to be a kernel point of X with respect to S. It is clear that s * is also a kernel point of S with respect to X. (x * , s * ) is said to be a pair of kernel points in X × S. The set of all and only kernel points of X is called the kernel set of X with respect to S. Similarly, the kernel set of S with respect to X contains only the kernel points of S. Figure 1 shows an example of 2D point sets S and X, the points reaching the minimum distances, the kernel points, the kernel sets, and the Voronoi cells.
Application to localization nanoscopy. RMSMD can be applied to localization nanoscopy imaging in various situations. S can be a set of emitter locations and X is a set of locations estimated from a sequence of data frames generated by the emitters. Given a set of emitter locations, RMSMD can evaluate the quality of a nanoscopy image estimated by a localization algorithm, and therefore evaluate the performance of the algorithm. The evaluation can be based on all the locations estimated from the entire sequence of acquired data frames. There are two categories of localization algorithms in the literature. In the first category, a localization algorithm like the Bayesian localization algorithm 10 jointly uses the entire sequence of data frames to estimate one location for each emitter, and therefore M and N are in the same magnitude. In contrast, belonging to the second category, most localization algorithms 1-3 localize the activated emitters frame by frame independently and do not identify an estimated location with a particular emitter. Since each emitter is activated a number of times in the sequence of data frames, these algorithms produce multiple estimated locations for each emitter. The nanoscopy image consists of the estimated locations from all data frames and therefore the number of estimated locations is much larger than the number of emitters, i.e., N ≫ M. For example, if each emitter is activated K times on average, then N ≅ KM on average where K can be as large as 50 for some kinds of emitters 11 . For the second category of localization algorithms, RMSMD can also be applied to evaluate the quality of emitter locations estimated from each single frame. In this case, S and X are the set of locations of the activated emitters in a frame and the set of their estimates, respectively, and then N and M differ slightly.
A sample drifting may cause additional translation and rotation for the entire set of estimated locations X and then the RMSMD D(X, S) might be unnecessarily increased. However, such an increase can be eliminated by performing the converse translation and rotation for the entire set of estimated locations. Because of this, throughout we assume that such additional translation and rotation have been removed from X.

RMSMD Properties
A number of properties of RMSMD are presented in this section. Some of them are common to a distance metric just as RMSMD, and some others are unique to RMSMD.
Properties common to a distance metric. We present the following properties of RMSMD that are common to a distance metric. The proof of them is straightforward and is omitted. Properties of kernel points. In addition to the properties common to a distance metric, RMSMD also has the following unique properties. Property 6 (A pair of kernel points): (x * , s * ) is a pair of kernel points if and only if they are located in each other's Voronoi cell, that is, x * ∈ V(s * ) and s * ∈ V(x * ). x 1 x 2 x 3 x 4 x 5 x 6 x (nm)  Property 7 (Removal of far non-kernel points): Suppose that X′ ⊂ X is a subset of non-kernel points such that (i) for any x ∈ X′, (ii) the average of square minimum distances − ∈ s x min s S 2 for all x ∈ X′ is larger than the average of the remaining square distances in D(X, S). Then − ′ < D X X S D X S ( , ) ( , ). If X * and S * are kernel sets for each other, then i i is the ith pair of kernel points. Property 8 (Kernel points): Let X * ⊆ X and S * ⊆ S be the sets of kernel points of X and S, respectively. If the average of square distances between each pair of kernel points (x * , s * ) ∈ X * × S * is smaller than the average of the remaining square distances in D(X, S), ( , ). Property 7 indicates that removing a subset of non-kernel points from X that are far away from the kernel points of S can decrease RMSMD. Property 8 further indicates that all non-kernel points can be removed to reduce RMSMD if the average of square distances between each pair of the kernel points is smaller than the average of the remaining square distances. Property 7 and Property 8 are proved in the Appendix.
For example, in Fig. 1 the MSMD between S and X is  It confirms Property 7 that removing X′ from X, we obtain The kernel sets are X * = {x 1 , x 4 , x 6 } and S * = {s 1 , s 2 , s 4 }. It confirms Property 8 that Looking at Fig. 1, it is clear that X − X′ is a better fitness of S than X, and X * and S * fit better to each other than X and S do, which are predicted by the corresponding RMSMDs. In return, it implies that RMSMD is a rational measurement of mutual fitness between two sets of points.
In localization nanoscopy, a majority of localization algorithms in the literature 1,3 localize emitters independently frame by frame and produce a number of estimated locations for each emitter. Most of the estimated locations are far away from the ground-truth locations and therefore can be removed to reduce RMSMD and improve the quality of nanoscopy image according to Property 7. As implied by Property 8, keeping only the kernel estimated locations can minimize RMSMD. Hence, an algorithm that without knowing the ground-truth locations can identify a number of non-kernel estimated locations with certain probability is desirable.

Accuracy, precision, recall, and JAC
Accuracy, precision, recall, and JAC are currently employed in the literature as quality metrics for localization nanoscopy images 1,3,8,9 . To demonstrate the advantages of RMSMD over these four metrics, in this section we investigate their ambiguity, discontinuity, and failure in distinguishing qualities of different images in certain conditions. Accuracy, precision, recall, and JAC all are defined on the basis of the true-positive regions of emitter locations s j 's. The area within a circle centered at s i is considered as the true-positive region of s i and the rest is considered as its false-positive region. The radius of the circle is equal to the FWHM of PSF in an optical system. An respectively, where a i for i = 1, …, |TP| is the distance from a true-positive point x k to a nearest s j . According to their definitions, accuracy, precision, recall and JAC have the following substantial drawbacks. The FP and FN ambiguities imply that these four metrics fail to distinguish the qualities of different nanoscopy images in certain conditions.
Remarks (Estimation vs. detection): In a statistical decision of the value of a variable through observations, there are two basic problems: one is estimation and the other is detection 12 . In an estimation problem, the variable with an unknown value is defined on a continuous support and so is continuous and can take uncountably infinitely many possible values. A device to make a decision on the unknown value in an estimation problem is called an estimator. On the other hand, in a detection problem, the variable with an unknown value can only take a finite number of possible values. A device to make a decision on the unknown value in a detection problem is called a detector. In a binary detection problem a detector can only choose one of two possible values, say true or false. In localization nanoscopy, since the emitter locations are defined on n  , i.e., a real space with dimension n = 2 or 3, the localization of an emitter is an estimation problem but not a detection problem. For the same reason, how well X (a set of estimated emitter locations) and S (the set of ground-truth emitter locations) mutually fit to each other is a problem to measure the performance of an estimator instead of a detector.
It is worthy to clarify what accuracy, precision, recall, and JAC are appropriate and inappropriate to measure. Precision and recall are suitable for measuring the performance of a detector in a binary detection problem. When applied to measuring the performance of emitter localization -an estimation problem, x i 's have to be divided into two groups, TP and FP. To this end, in refs 1,3 the circle centered at s i with the radius of PSF FWHM in an optical system is used to divide the space into the true-positive and false-positive regions. Such a division depends on the optical system in an experiment, and so is somewhat subjective and arbitrary, and is not universal as well. JAC is originally suitable for measuring the fitness of two continuous sets consisting of uncountably infinitely many points. However, X and S in localization nanoscopy are two discrete sets each consisting of countably many locations, and so their mutual fitness is inappropriate to be measured by JAC. The original accuracy is suitable to measure the global error between an estimated value and its ground-truth value. When accuracy is applied to localization nanoscopy, the partition set X i of X for each ground-truth location s i must be known. In practice, however, only the entire set X of estimated emitter locations is known but its partition sets X i 's are unknown. To circumvent this problem, in refs 1,3 those x j 's located in the true-positive region of s i are considered as estimates of s i and are taken into account in calculation of accuracy in Eq. (8); and the other estimated locations that are not located in the true-positive region of any s i are ignored. Because of this, the accuracy based on the true-positive region cannot properly measure the overall fitness between X and S. Even when the partition X i 's were known, the accuracy in Eq. (8) would only measure the error between a ground-truth location and its estimates and cannot properly measure the average, local, and mutual fitness between X and S that truly needs to be measured.
In summary, precision and recall are not suitable for performance evaluation of an estimator. The modified JAC and the modified accuracy still fail to properly measure the mutual fitness between two discrete sets of points. In contrast, as a universal and objective metric, RMSMD can properly measure the average, local, and mutual fitness between two sets of points.

Numerical Examples
In this section, the performances of RMSMD and the four metrics of accuracy, precision, recall, and JAC are compared by two numerical examples.  Fig. 2(b)-(d).
As shown in Fig. 2(b), RMSMD continuously properly measures the fitness between X and S in all cases as X continuously changes.
In contrast, as shown in Fig. 2(c) and (d), accuracy, precision, recall and JAC present drawbacks. First, they are discontinuous with respect to d. Second, when x 8 , …, x 14 are in the false-positive region, accuracy is undefined, and precision, recall, and JAC are constantly equal to zero. That is, they all fail to distinguish the qualities of continuously changing X as x 8 , …, x 14 move along the corresponding segment of trajectories. Third, when x 8 , …, x 14 are in the true-positive region, precision, recall, and JAC are always equal to constants 50%, 100%, and 50%, respectively, and fail to distinguish the qualities of continuously changing X. Finally, in the particular case when x 8 , …, x 14 are coincided with s 1 , …, s 7 , respectively, accuracy is equal to zero, which implies an illusive perfect fitness between X and S but actually x 1 , …, x 7 are not taken into account. These drawbacks are inherently incurred by the FP and FN ambiguities, FP and FN discontinuities, and inappropriateness of accuracy, precision, recall and JAC as quality metrics. x 4 x 5 x 6 x 7 x 8 x 9 x 10 x 11 Nanoscopy imaging with an unbiased Gaussian information-achieving estimator. We simulate an experiment of 2D localization nanoscopy imaging where an unbiased Gaussian estimator that achieves the Fisher information and CRLB localizes activated emitters frame by frame independently. The RMSMD of an information-achieving nanoscopy image can benchmark localization nanoscopy images obtained by other algorithms. The data frame model, Fisher information matrix, and CRLB are defined and presented in detail in ref. 7 to which the reader is referred. Fig. 3(a), M = 250 emitters are located on a 2D helix with the adjacent emitter distance equal to 25 nm. Each data frame is generated according to the model in ref. 7 . The optical system has a numerical aperture of n a = 1.4 and the fluorescence wavelength of emitters is λ = 540 nm. Then the standard deviation of an Airy PSF can be approximately calculated by Eq. (49) in ref. 7 and is equal to σ = 78.26 nm, and accordingly σ = = . ln FWHM 2 2 2 184 28 nm. A Gaussian PSF with zero mean is used to approximate the Airy PSF with the same standard deviation σ. An activated emitter in a frame emits a number of photons that is Poisson distributed with the mean of I photons per second. Each photon is detected by the camera with a probability of η. The mean number of photons that are emitted from an activated emitter and are successfully detected by a camera is equal to Iη = 300000 photons per second. The frame rate is 1/Δt = 100 frames per second. Then the mean number of detected photons is equal to IηΔt = 3000 per frame per emitter. The additive Poisson noise and Gaussian noise are presented in the data frames 7 . The signal to Poisson noise ratio (SPNR) and the signal to Gaussian noise ratio (SGNR) are SPNR = 0.2 and SGNR = 0.3 µm 2 per emitter, respectively. The frame size is L x × L y = 2048 × 2048 nm 2 and the pixel size is Δ × Δ = × 128 128

Data frame model. As shown in
x y nm 2 , and then each frame has K x × K y = 16 × 16 pixels. The photon emissions of all emitters, Poisson noise, and Gaussian noise are all independent. Emitter activation process. The activation process of an emitter can be modeled by a Markov chain 10,13 . We consider a more general chain than those of refs 10,13 . Specifically, each emitter in a data frame is independently activated and the state of activation follows a Markov chain taking states 0, 1, 2, 3, 4. State 0 means no activation and the other four states mean activation. The transition probability from state i to state j is denoted by r ji . In the simulation, > r 0 00 is variable, and r 01 = 0.5, r 02 = 0.7, r 03 = 0.8, and r 04 = 1 are fixed; correspondingly, r 10 = 1 − r 00 , r 21 = 1 − r 01 , r 32 = 1 − r 02 , and r 43 = 1 − r 03 , and all other transition probabilities are zero. It is clear that the Markov chain is irreducible, aperiodic, and positive recurrent and therefore has a unique stationary probability distribution 14  To cover a broad range of emitter density, the average density of activated emitters in a data frame is purposely set up higher than that in practice 1,3 .
Unbiased Gaussian estimator achieving the Fisher information. In each data frame, a random number K of emitters are activated as shown in Fig. 3 2 are called the CRLB 6,7 that are the smallest variances possibly achievable by all unbiased estimators for x 1 ,y 1 ,…,x K ,y K , respectively. An unbiased estimator that achieves the Fisher information, i.e., whose Fisher information matrix is F, also achieves the CRLB.
To generate an unbiased Gaussian estimator θˆ achieving the Fisher information, we consider a Gaussian estimator θˆ with mean θ and covariance matrix F −1 , i.e. θ θ −~N F ( , ) 1 . It is straightforward to verify that the Fisher information matrix of θˆ is equal to F. Therefore, θˆ is unbiased and Gaussian and achieves the Fisher information and the CRLB as well. Denote by F −1 = UΛU T the eigendecomposition of F −1 where Λ is the diagonal matrix of eigenvalues and U is the matrix of corresponding eigenvectors of F −1 . Let g be a 2K-dimensional standard Gaussian random vector with zero mean and unit covariance matrix. Then an unbiased Gaussian estimator of θ that achieves the Fisher information and CRLB is obtained by It is clear that θ θ −~N F ( , ) 1 . Accordingly, the ith 2D vector of θ θ θ = .ˆ( ) , , i.e., σ xi 2 and σ yi 2 , are the CRLBs of x i and y i , respectively. Therefore, θˆi is an unbiased Gaussian estimator achieving the CRLB of θ i . A realization of θˆ can be obtained with given a realization of g, say a pseudorandom vector generated by the embedded function randn in MATLAB. The realization of θˆi is one of the estimated locations from a data frame and is shown in the final nanoscopy image as presented in Fig. 4.

95%-probability region.
It is interesting to know where θˆi is mostly located. To this end, we define the 95%-probability region of θˆi as the region where θˆi is mostly located and the total probability is 0.95. Since θˆi is Gaussian distributed with mean θ i and covariance matrix − F i 1 , it is clear that the 95%-probability region is given by the ellipsoidˆθ i is Gaussian distributed with zero mean and unit covariance matrix, v i is Chi distributed and its number of degrees of freedom is equal to n where n is the dimension of v i . Therefore, P(n/2, R 2 /2) = 0.95 where P is the regularized gamma function, the cumulative distribution function of v i . When n = 2, P becomes the Rayleigh distribution 15  , and can be drawn in a two-dimensional image. The larger the 95%-probability region, the larger the estimation error of θˆi.
For the activated emitters in Fig. 3(b), their 95%-probability regions of the information-achieving estimators are shown. As can be seen, the emitters located close to each other with their PSFs severely overlapped are suffered from the significant inter-emitter interference as discussed in ref. 7 and have a large 95%-probability region, and their estimation errors in x and y directions may be considerably correlated due to a large difference between the two eigenvalues. A data frame generated according to the data frame model is shown in Fig. 3(c). By means of Eq. (12) the estimated emitter locations by the unbiased estimator achieving the Fisher information are shown in Fig. 3(d). The nanoscopy image consists of all the emitter locations estimated from the N data frames as shown in Fig. 4(a)-(d). Since all activated emitters are localized frame by frame independently, the nanoscopy images achieve the Fisher information and CRLB in the sense that the estimated emitter locations achieve the Fisher information and CRLB in each corresponding data frame.
Indexed by the mean number of activated emitters per frame M a , the RMSMD, accuracy, precision, recall, and JAC of the nanoscopy images are shown in Fig. 5(a) and (b). As the mean number of activated emitters per frame increases, the RMSMD increases exponentially fast, which confirms the result that the CRLB exponentially increases as the emitter density increases 7 , indicating that the quality of the nanoscopy image decreases fast. The fast increase of both RMSMD and CRLB is due to the estimated locations spread away from the ground-truth locations caused by the increasing severity of the inter-emitter interference. Subjectively, the visual quality of the nanoscopy image indeed deteriorates fast accordingly. In contrast, the accuracy, precision, recall, and JAC all change slightly and are insensitive to the quality change. The result demonstrates the superior sensitivity of RMSMD to a quality change over the other four metrics.

Conclusions
We propose RMSMD as a quality metric for a localization nanoscopy image that consists of estimated emitter locations in comparison with their true locations. RMSMD depends only on the two sets of points regardless of how they are obtained. RMSMD measures how well two sets of points averagely, locally, and mutually fit to each other. As a distance metric, RMSMD possesses the properties common to distance metrics as well as its own unique properties. The four metrics of accuracy, precision, recall, and JAC inherently present ambiguity, discontinuity, and inappropriateness and fail to distinguish the qualities of different localization nanoscopy images in certain conditions. In contrast, RMSMD presents the advantages of universality, objectiveness, continuity, and sensitivity to a quality change over these four metrics. The unbiased Gaussian information-achieving estimator is a benchmark for practical localization algorithms. The RMSMD of an information-achieving localization nanoscopy image in addition to its visual quality benchmarks the quality of localization nanoscopy images. As a universal and objective metric, RMSMD can be broadly employed to evaluate the quality of localization nanoscopy images and the performance of localization algorithms and in general to measure the mutual fitness of two sets of points as well.
Code availability. The customer code of Matlab Version 7.5.0.342 (R2007b) that produces the figures is available online as supplementary material.