Image Anal Stereol 2000;19:209-214 Original Research Paper ESTIMATING NUMBER DENSITY Nv - A COMPARISON OF AN IMPROVED SALTYKOV ESTIMATOR AND THE DISECTOR METHOD Ashot Davtian1, Ute Hahn2, Joachim Ohser3, Dietrich Stoyan1 Institut fu¨r Stochastik, TU Bergakademie Freiberg, 09596 Freiberg, Germany; 2Department of Mathematics and Statistics, The University of Western Australia, Nedlands WA 6907, Australia; 3 Institut fu¨r Techno- und Wirtschaftsmathematik, 67663 Kaiserslautern, Germany (AcceptedMarch 17, 2000) ABSTRACT Two methods for the estimation of number per unit volume NV of spherical particles are discussed: the (physical) disector (Sterio, 1984) and Saltykov’s estimator (Saltykov, 1950; Fullman, 1953). A modification of Saltykov’s estimator is proposed which reduces the variance. Formulae for bias and variance are given for both disector and improved Saltykov estimator for the case of randomly positioned particles. They enable the comparison of the two estimators with respect to their precision in terms of mean squared error. Keywords: disector, number density, planar sections, precision, Saltykov estimator INTRODUCTION Estimating the mean number of particles per unit volume NV is a fundamental stereological problem. There are two main approaches, the choice of which is mostly driven by practical considerations. In cases where the structure of interest can be investigated as a three-dimensional structure, as e.g. in confocal microscopy, the (optical) disector is the most convenient method, and it gives unbiased estimates of NV (Sterio, 1984; Howard etal, 1985). Therefore it is very well established in biological sciences. However if the matter is opaque, the physical disector has to be used which is based on parallel sections a known distance apart. Then it is necessary to identify particles that are hit by both sections, which requires perfect registration of the relative position of the two images. Moreover a bias may be introduced because particles that are smaller than the distance between the two section planes cannot be detected. Both optical and physical disector imply a lot of manual work of qualified assistants. Often it is practically very difficult or even impossible to make parallel sections, as of hard and brittle material. This is a reason why classical stereological methods using single sections are traditionally preferred e.g. by materials scientists. These methods are based on assumptions about the shape of the particles. Perhaps the most popular estimator is Saltykov’s method (Saltykov, 1950; Fullman, 1953). It assumes that the particles are spherical, an assumption which in engineering is very often quite realistic. Saltykov’s estimator requires only to measure the diameters of the (circular) particle profiles. In many cases this can be accomplished by automated image analysis. The estimator is often blamed for its theoretically infinite variance, which is caused by very small intersection circles. As shown in this paper, it can easily be improved by binning these small diameters. Then its variance is finite at the cost of a negligible bias. Up to now, hardly anything has been known about the statistical properties of the physical disector and the improved Saltykov estimator. It is clear that the variance of the estimators largely depends on the degree of regularity of the investigated structure, however the relationship is very complex, and in practical applications it is difficult to investigate the spatial variability of particle arrangement. We give formulae for variance and bias of both methods based on the assumptions that the particles are spherical and completely randomly positioned. Accuracy of the estimators is computed for various distributions of sphere diameters. The results lead to practical recommendations for the observation window size which is necessary to obtain a given mean squared error. ESTIMATORS OF Nv A well-known classical planar-section estimator N ˆ V S ofNV (Saltykov, 1950; Fullman, 1953) is based on the observation of section circles with centres in a planar window W of area A(W). Their (random) number is n, their diameters are s1,...,sn. Then 209 Davtian A et al: Estimating NV - improved Saltykov estimator vs. disector This estimator is unbiased, given the section plane is randomly positioned or, in terms of the model-approach, the system of spheres is homogeneous, i.e. can be described by a stationary germ grain model (Stoyan et al, 1995). Unfortunately, the variance of NV is infinite due to the occurence of very small section circles. In practice microscopical image resolution sets a natural lower limit to the diameters that can be observed, i.e. there is a lower bound for the si. But still the variance of NV S will be very large. Thus it is natural to modify the estimator such that very small section circles are not measured but only counted. To this purpose a lower diameter e is introduced such that all section circle diameters larger than e are measured while all smaller ones are set equal e/2. With the notation s*{e) = Ui2 if si i I s, otherwise the modified estimator is N S* 2 n 1 nAi 1s*(e (2) Clearly, NV S must have a bias which increases with increasing e. As shown below, the bias is rather small and the variance of NV S is finite, decreasing with increasing e. The physical disector estimator NV (Sterio, 1984) is based on counts made on two parallel section planes of distance t. The number Q~ of particles which appear in a counting frame in the reference plane but do not hit the look up plane has to be determined. With t denoting the distance of the two planes and A(W) the area of the counting frame, N ˆ D Q- tAWy (3) This estimator is unbiased only if there are not any particles with diameter smaller than t. Otherwise such small particles may be situated between the two planes and therefore are overlooked, which results in a negative bias. Furthermore, the disector requires the distance t to be precisely known and that it is possible to decide whether profiles in the two planes belong to the same particle or not. The same assumptions about homogeneity of the particle system or randomization of the section planes are made as for the estimators NV andNV S. STATISTICAL PROPERTIES OF THE Nv ESTIMATORS Under the homogeneity assumptions made in the previous section, the means and biases of NV S and NV can be calculated quite easily. In this section only results are presented, for sketches of their derivation see the Appendix. Means (and biases) of the estimators depend on the distribution function DV of the sphere diameters. In the Appendix it is shown that EN ˆ V S =NV 4 ne 2 ("" e — / arccos(-)DV(d d) K Je d du (4) and enV D=nv[1 t J0 (t-d)DV(d d) (5) Using these formulae the relative bias B' = (ENV—NV)/NV was calculated for the two estimators. Five diameter distributions were considered, all with mean one: constant diameters, uniform (on [0.7,1.3]) and triangular (on [0.1,1.9]) as examples of bounded distributions, as well as Rayleigh distribution and lognormal distribution (with variance 0.1). The corresponding density functions are depicted in Fig. 1. 1.8 .2 1.6 | 1.4 I 1.2 S-H H 1 I 0.8 -o 0.6 0.4 0.2 0 o 'm C Q i/ * i' 1 "C : /1 '¦> i R// 'i .• / ¦ .' -4/1 '¦ i v i / •' ' / •' ' / .' / 1 1 1 1 i \ X \ X i \ v\ 's '. \. 1 "*'-. ^""""----^ 0 0.5 2.5 1 1.5 2 Sphere diameter d Fig. 1. Density functions of the diameter distributions. U - uniform, T - triangular, R - Rayleigh, L -lognormal. For the parameters see text. ( C - line representing constant diameters = 1.) Since the relative bias is dimensionless, the results also hold for diameter ¯ distributions with the same shape but with mean dv different from one. Tothis end ¯ B' is regarded as a function of the ratio e : dV or t: dV, resp. 210 Image Anal Stereol 2000;19:209-214 a) > VarN V 0.4 0.6 Threshold e Q (b) I 0.3 0.2 0.1 -0.1 P4 -0.2 -0.3 -0.4 0.2 0.4 0.6 0.8 1 Section distance t 1.2 1.4 Fig. 2. Relative bias (a) of the improved Saltykov estimator and (b) of the_disector as a function of the ratios e : dV or t : dV, resp. for the diameter distributions in Fig. 1. dV G [0,1] and (b) BD Fig.2 shows (a) BS for e fort:dVe [0,1.4]. The absolute value of the relative bias B S of the improved Saltykov estimator is very small (between 0.0004 and 0.0026) for all considered distributions if £ < 0.2. It still does not exceed 0.0082 for e < 0.4. Therefore a ratio "e : mean diameter” of 0.2 or even 0.4 can be recommended for practical use. The disector is unbiased if the distance t is smaller than the smallest sphere diameter, e.g. for t < 1 in case of constant diameters. In general the bias of the disector depends very much on the form of the diameter distribution. In contrast to mean or bias, the variance of the estimators also depends on the variability in spatial arrangement of the particles. In general the relations are very complex. However it is possible to give formulae for the case of independently random positioned sphere centres and independent diameters, i.e. for the Poisson case. The variance of NV S is given by NV A(W) k2 h d Ve 1)DV(d d)+ (6) 16 K2£2 J0 (DV(\/u2 + £2) -DV(u))du and the variance of N D is VarND Nv A{W) see the Appendix. V .1(i-1 t 1 t (t-d)DV(d d)), (7) In parallel with relative bias, relative variance V = VarNV ¦ (dVA(W))/NV was calculated for both estimators as a dimensionless quantity. Fig. 3(a) shows VS' for e G [0,1]. The variance does not very much depend on the distribution of the diameters. As to be expected, it decreases with increasing e. For e = 0.2, VS* takes values between 1.75 (constant diameters) and 2.18 (Rayleigh distribution), and for e = 0.4, VS' G [1.48,1.71]. The relative variance VD of the disector method is shown in Fig. 3(b). Obviously it depends even less on the type of diameter distribution than the variance of NV S does. It is quite large (VD > 3.25) for t < 0.3, but it decreases with increasing t because then the volume of the three-dimensional counting box increases which is used for the disector estimator, i.e. the sample size increases. PRECISION OF THE Nv ESTIMATORS Performance of an estimator is usually expressed by the mean squared error (MSE), that is the average squared distance of the estimate to the true estimated parameter. MSE is a combination of variance and bias: MSE' = E (NV-NV)2 = VarNV + (biasN V)2. For any number density NV and any window size A(W), the MSE of the NV estimators can easily be calculated from relative variance V and relative bias B': MSE' Nv -V+NV 2(B)2=NV 2(V_ + (B)2), (8) dVA(W) where En = dVNVA(W) is the mean number of section circles counted in the observation window or counting frame, respectively. This enables the comparison of the precision of the estimators NV S and NV using the formulae and diagrams of the previous section. Clearly, MSE depends on the diameter distribution and on NV 211 Davtian A et al: Estimating NV - improved Saltykov estimator vs. disector (a) °L > is Pi (b) Pi 0.4 0.6 Threshold e 3.5 \ 3 \ 2.5 V 2 \v 1.5 1 0.5 C^_ 0 0.2 0.4 0.6 0.8 1 Section distance t 1.2 1.4 Fig._3. Relative variance (a) of the improved Saltykov estimator and (b) of the disector as a function of the ratios e : dV or t: dV, resp. for the diameter distributions in Fig. 1. itself, two factors that are intrinsic to the investigated structure and cannot be altered. Yet the experimenter can influence MSE by the choice of observation window areaA W) and design parameters e or t, resp. A popular rule of thumb advises to chose the sample size — here represented by the window area — such that the variance part of MSE is not smaller than the bias part. This means, in terms of (8), that V /En should be greater or equal to (B')2. The two parts are of equal size if the observation window or counting frame size is chosen to beA(W)opt = V/(dVNV(B')2). As an example, two cases are studied of a structure with independent randomly positioned spheres of constant and independent Rayleigh distributed diameters, resp. The number density is N V 1000 mm -31 and the mean diameter is dV = 0.02 mm. For the disector, the section distances t = 0.02 and t = 0.01 [mm], and for the improved Saltykov estimator, e = 0.008 and e = 0^004 [mm] are considered, corresponding to ratios t : dV of 1 and 0.5 resp., and e : dV of 0.4 and 0.2, resp. Table 1 shows relative bias and variances, as well as the MSE for a window of area 10 (“MSE10”), and the area A5000 of a window such that MSE = 5000. The absolute bias of the improved Saltykov estimator is very small for any of the considered distributions. Its contribution to MSE is practically negligible. For the physical disector, the bias strongly depends on the diameter distribution. If the structure has a relatively high proportion of small particles, a small bias can only be obtained with a small section distance t. But then the variance is very large, again leading to a large MSE. E.g. in the case of the Rayleigh distribution, NV with t = 0.5dV has both a larger relative bias_and a larger relative variance than NV S with e = 04dV, see Table 1. As the bias of NV S (0.4) is negligible compared to the variance and as variance of NV increases with decreasing section distance, NV Si04) is better than any disector with distance t < 0.5dV. Generally spoken, the improved Saltykov estimator is particularly useful if the diameter distribution has a large range. Taking into account that the disector uses two section planes, the total window area which is required Table 1. MSE10 = MSE [mm~6] for a window of area 10, and A5000 = area [mm2] of a window such that MSE =5000, of disector and improved Saltykov estimator for different parameters t and e, resp. NV= 1000 mm-3], dV = 0.02 mm. For further details see text above. constant diameters Rayleigh distribution B" V MSE10 A 5000 B' V MSE10 A 5000 NVD: t d_y= 1 0 1.000 5000.0 10.0 -0.2101 0.7899 48087.9 * t dV = 0.5 0 2.000 10000.0 20.0 -0.0618 1.8765 13197.6 79.2 NS': e : dy = 0.4 0.0038 1.4809 7418.5 14.9 -0.0080 1.7118 8622.3 17.3 e : dV = 0.2 0.0004 1.7479 8739.8 17.5 -0.0010 2.1811 10906.7 21.8 A5000 does not exist, MSE > 5000 for any window area. 212 Image Anal Stereol 2000;19:209-214 to achieve the same precision is smaller for NV S even in the case of constant diameters. However measuring diameters in the case of NV S means acquiring more information of the window than bare counting as in the case of NV, and therefore perhaps a little higher effort to evaluate the sections. A preliminary report was presented at the Xth International Congress for Stereology, Melbourne, Australia, 1-4 November 1999. REFERENCES FullmanRL (1953). Measurement of particle size in opaque bodies. J Metals 5:447-52. Howard CV, Reid S, Baddeley A, Boyde A (1985). Unbiased estimation of particle density in the tandem scanning reflected light microscope. J Microsc 138:203-12. Saltykov SA (1950). Introduction to Stereometric Metallography (Russ.) Erevan: Publ Acad Sci Armen SSR. Sterio DC (1984). The unbiased estimation of numbers and sizes of particles using the disector. J Microsc (1984) 134:127-36. Stoyan D, Kendall WS, Mecke J (1995). Stochastic geometry and its applications. Chichester: John Wiley & Sons. APPENDIX: DERIVATION OF ESTIMATOR MEANS AND VARIANCES Following the notation in Stoyan et al. (1995, p. 355), the particle system is represented by a spatial germ grain model (or marked point process) *¥ V = {[xn,yn,zn;dn]} of spherical particles, where (xn,yn,zn) is the centre of the nth sphere and dn is its diameter. The corresponding unmarked point process of particle centres has density NV. In order to calculate the mean of NV, we translate the number Q~ into that notation. A sphere of diameter d with centre in (x,y,z) intersects both reference and look up plane if z— d 2 < 0 < z+d 2 < t. Then it produces a section circle of radius r = ^(d/2)2 — z2. The circle is counted in the counting frame if its centre falls inside the frame shifted by the vector (r, r). Hence N D Q- tA(W) 1 W > [x,y,z;d\eWV (x,y)- 1( z--<0 [x,y,z;d\eWV ¦1(-d/2 e)Vd2-4z2+ + 1(Vd2-4z2e) 1(Vd2-4z2