Image Anal Stereol 2012;31:99-108 Original Research Paper LOCAL STEREOLOGY OF EXTREMES Zbynek Pawlas^3 Department of Probability and Mathematical Statistics, Faculty of Mathematics and Physics, Charles University, Sokolovska 83, 18675 Prague, Czech Republic e-mail: pawlas@karlin.mff.cuni.cz (Received November 28, 2011; revised April 2, 2012; accepted April 4, 2012) ABSTRACT Local stereology uses information obtained from central sections passing through a reference point of the particle. The aim of this paper is to investigate the prediction of extremes of shape and size parameters based on the central sections. We consider the particle population formed by spheroids (either prolate or oblate) and assume that the reference point is the centre of the spheroid. a relation between shape and size parameters of the particles and their planar sections is derived and consequently stability properties of the domain of attractions are proved. Keywords: extremes, local stereology, maximum domain of attraction, shape and size parameters, spheroids. INTRODUCTION Local stereological methods form quite modern branch of stereology, see Vedel lensen (1998) for a comprehensive exposition. They require that we associate a reference point to each particle and accomplish sectioning through this reference point. We speak about central sections or local probes. Information about the particle population is then extracted from these local probes. In contrast to traditional sampling, where the system of particles is sectioned by an isotropic uniform probe, the local sample is representative for the particle population. This fact was exploited in Pawlas et a I. (2009) for the development of a statistical procedure for obtaining information about particle size distribution from central sections without specific assumptions about particle shape. The motivation for local stereology comes from the study of biological tissues. The particles are cells and the centres of cell nuclei or nucleoli serve as reference points. Local techniques are most easily implemented if optical sectioning is available. One of the possible applications in stereology is to estimate extremes of particle parameters from the observation of test probes of lower dimension. This field is referred to as stereology of extremes (Benes and Rataj, 2004). In practical applications it is often important to analyze extremes of the particle parameters. For example, the damage of materials is related rather to extremal than mean characteristics of microstructure. The application of stereology of extremes to metallurgy is discussed e.g., in Takahashi and Sibuya (2002) or Bortot etal, (2007). So far in the literature concerning stereology of extremes, only isotropic uniform random probes were considered. In Wicksell's corpuscle problem (Wicksell, 1925), the prediction of maximum size of spherical particles was studied by Drees and Reiss (1992) and Takahashi and Sibuya (1996; 1998; 2001) while the behaviour of minimum size was investigated in Kotzer and Molchanov (2006). An extension to spheroidal particles leads to the prediction of not only extremal size but also extremal shape. This was considered in Hlubinka (2003a;b; 2006a); Hlubinka and Kotz (2010) for oblate spheroidal particles and in Hlubinka (2006b ;c) for prolate spheroidal particles. The novel contribution of this paper is the study of stereology of extremes for spheroidal particles under different stereological sampling design. Our aim is to investigate the stereological estimation of the tail of particle shape and size distribution based on the local probes. We derive the relation between the particle parameters and the parameters of particle central section (also called profile). Subsequently, we show how the limiting distribution of the extremal particle parameters is related to the limiting distribution of the extremal profile parameters. It turns out that the distribution of profile parameters belongs to the same domain of attraction as the distribution of particle parameters. Particular attention is devoted to the population of spheroids (either prolate or oblate). If we consider a spheroid with two semiaxes of length a and one semiaxis of length c, then its planar section through the centre is an ellipse with semiaxes of lengths a and d, where inin(Y/.r) < d < inax(Y/.r). In other words, the length of one semiaxis is directly observed at the central section while the second semiaxis appears shorter (for prolate spheroid) or longer (for oblate spheroid). Therefore, profiles have more circular shape than original particles projected to the plane perpendicular to one of semiaxes of length a. This paper is organized as follows. At first, we derive the relation between particle and profile shape and size parameters in the population of spheroidal particles. Then, we summarize basic facts from extreme value theory which will be needed to obtain main results stated in the subsequent section. We make the inference about the extremal domains of attraction under local sampling design. We conclude with an example showing the behaviour of maximal shape parameter in the population of simulated oblate spheroids. PARTICLES Local stereological estimates are based on information collected on a section plane through a fixed reference point of the particle. In the present paper, we will concentrate on three-dimensional particles and two-dimensional section planes. By a particle K we understand a non-empty compact subset of R3. We associate a reference point a* to the particle K. A section plane through the reference point takes the form x + L, where L is a fixed two-dimensional linear subspace of R3. The information about the particle is deduced from the planar central section Kn (x + L). In this paper we pursue the modelbased approach. However, all results remain valid in the same way for the design-based approach. In the latter, the particle K is deterministic while the section plane is random. The class of all particles equipped with the Hausdorff metric forms a separable metric space. Thus, it is possible to define random particles, their distribution and independence. For more details on theory of random sets we refer to Molchanov (2005). We consider a population Ei,...,E„ of random particles with reference points X\,..., Xn. The particles are not observed directly, only their planar profiles Ei n (Xi + L) are available. An illustration for two-dimensional particles and one-dimensional central sections is shown in Fig. 1. We assume that the particles Ei,...,E„ are independent and identically distributed. Denote by So a random particle with the same distribution as the Ei. We will refer to So as a typical particle. Let p be a real measurable function defined on the space of all particles. It will be used to describe shape or size parameters of particles. Our aim is to investigate the Fig. 1. A two-dimensional illustration of local sectioning applied to the particle population under study. tail behaviour of the typical particle parameter p(So). Since the parameters p(Zi),...,p(Zn) describing random shape or size of particles cannot be directly observed, we have to use the particle profile parameters q{Si n (Xi + L)),...,q{E„ n [Xn + L)) derived from local probes. In what follows, we will consider uniform randomly oriented (isotropic) particles. If the particles cannot be regarded as isotropic, we may instead randomize the orientation of the section plane (use isotropic subspace L) and work under design-based setting. If So is a triaxial ellipsoid with semiaxes of lengths A, B and C (A > B > C) and the corresponding reference point Xq is the centre of the ellipsoid, then the observed planar profile So n (Xo + L) is an ellipse. In particular, we will consider the case of two equal semiaxes, i.e., So is a spheroid. It can be either prolate (B = C, two equal minor semiaxes) or oblate (A = B, two equal major semiaxes). A possible choice for p describing shape of the typical particle is given by the shape factor This parameter is commonly used for spheroids, see, e.g., Benes and Rataj (2004); Hlubinka (2003a; 2006a;b). For spherical particles (A = C) we get 5 = 0, more elongated spheroids lead to larger value of S. In similar way, we can define the function q to quantify the shape of the planar profile. We will denote the profile shape factor by T. Each spheroid is uniquely determined by the centre X (reference point), the length A of major semiaxis, the length C of minor semiaxis, and the angles 0 and which describe its orientation. Let 0 be chosen as the angle between the norm vector of L and the major (for prolate case) or minor (for oblate case) axis, i.e., 0 is the latitude and is the longitude. In what follows we assume that the random vectors (A,C) and (0,) are independent. The profile observed on planar central section is an ellipse and it can be easily shown that it has semiaxes of lengths D A2C2 , C2sin20+ A2cos2© and C in the case of prolate spheroid and (1) D A2C2 A2 sin 0 + C2 cos2 0 and A for oblate spheroid. Shape factor of the ellipse obtained from central section is D2 Ssin20 C2" ~ 1 +5cos20 (2) for prolate spheroid and A2 r = —y — 1 = Ssin2 0 for oblate spheroid. Note that T < S in both cases. It means that the profile shape factor is always smaller or equal to the particle shape factor. Equality occurs if and only if the section plane is perpendicular to the plane spanned by two semiaxes of the same length. Since the spheroid So is assumed to have isotropic orientation, the density of 0 is /©(0) = sin0, 0 t Ssin20 1+Scos20 rn/2 Jarcsin y/í(l+s)/.s(l+í) i r ls-t a/T+7 Jt v s for prolate case and 1 —Fr(t) =P(Ssin20 > t) r<*> i-arcsin\ftfs smdddFsids) dFs(s) , t > 0 , sindddFsids) s-t dFs(s) , t > 0 , for oblate case, where Fs is the distribution function of S. It will turn out to be useful to rewrite these formulas using integration by parts. For prolate spheroid it yields 1 ~FT(t) = 1 Fs(s) s/T+t y/l+tJt 2 S3/2s/^t d.s which lor i > 0 becomes 1 ~FT(t) = t f00 1 -Fs(s) VTTtJt 2s3/2V^t ds, t> 0. (3) For oblate spheroid we end up with the following formula 1 -FT(t)=t r i -Fs(s) 't 2 s3/2V^t ds, t> 0. (4) We deal with the stereological unfolding problem. The estimation of unknown particle shape distribution function Fs based on the estimator of Fj is an ill-posed problem. Based on the extreme value theory we are going to study the tail behaviour of shape factor. It is natural to use the semiaxes lengths as the size parameters. Since one semiaxis can be recovered from central section, we condition on the knowledge of its length and we are interested in the other semiaxis length, that is, A for prolate spheroids, C for oblate spheroids and D for profiles. Similarly as for the shape factors, we can derive the relation between the conditional distributions of particle size parameter and profile size parameter D. For prolate spheroids, we have from Eq. 1, l-FDlc(d I c) = P(D >d\C = c) A2(d2-c2) d2(A2-c2) sin 0 > /q2(rf2_c2 rf2(a2_c2 sin0d0,ÍA|C(d£71 c) a2 -d2 ~2-2 FA|c(dfl \c) , d>c. We apply the integration by parts and obtain 1 -FD\c(d | c) c(d2-c2) r a{l-FA\c{a\c)) : da , d Jd d > c . (5) Clearly, FD\c(d \ c) = 0 if d < c. For oblate spheroids we have FD\A(d I a) = P(D d2(a2 —C2) d pit/2 0 Jarcsin (a1-d1)c2 sin0d0Fcu(dc | a) d2(a fd a d2-c2 dV a2 — c2 Fc\A(dc I a) , d < a If Fcu(0 a) 0, then the integration by parts implies The normalizing constants can be chosen such that FD\A(d | a) a(a2 —d2) fd cFc{A(c\a) -de, d Jo ^^¿{a2 -c2)3/2 d < a . (6) Clearly, FDu(d \ a) = 1 if d > a. EXTREME VALUE THEORY Let Zi,..., Z„ be independent identically distributed random variables with common distribution function F. We are interested in the behaviour of the sample maximum Mn = max(Zi,... ,Z„). We say that F belongs to the maximum domain of attraction of a distribution function G if there exist normalizing constants cn > 0 and dn 0 and d'n such that lim — n^oo Cn i. iim dLJi = o. n^oo c then Eq. 7 holds with cn replaced by c'n and dn replaced by d'n. It is well-known that there are three possible non-degenerate limit distributions (Fisher-Tippett theorem, cf. Embrechts etal, (1997), Theorem 3.2.3), G belongs to the type of one of the following distribution functions (y > 0): - Frechet: Gi,r(z) = exp{-z~r}, z > 0, - Weibull: G2,r(z) = exp{-(-z)7}, z < 0, - Gumbel: G3(z) = exp{-e_=}, z € R. We summarize several results concerning the characterization of MDA(G), for more details we refer to de Haan (1975) or Embrechts et al, (1997), Section 3.3. Let coF = sup{z £ R : F(z) < 1} be the right endpoint of F and let F (u) = inf{z € R : F(z) > u}, ii G (0,1), be the quantile function. The distribution function F belongs to the maximum domain of attraction of G|.r if and only i!' oj/ oo and lim 1 -Fjiiz) 1 -F{u) z~r, for all z > 0 . cn = F (1 — l/n) and dn = 0. The distribution function F belongs to the maximum domain of attraction of Gi.y if and only if cop < oo and 1 -F{cof-uz) y hrn -——-- = zr, 1 — r ((Op — ii) for all z > 0 . The normalizing constants can be chosen such that cn = 0)/ —F (1 — l/n) and dn = oj/. The distribution function F belongs to the maximum domain of attraction of G3 if and only if there exists some positive function b such that lim It^COp- 1 -F{u + zb{u)) l-F(u) The auxiliary function to be differentiable on liin„ ,(0/ l/(ii) = 0, lim„ for all z G R • (8) b may be chosen (—00 ,coF) such that b{u)/u = 0 if (Dp = 00 and liin„ ,0)/ b(ii)/(o.)i- — u) = 0 if cop < 00. The normalizing constants can be chosen such that dn = — l/n) and cn = b{dn). Analogous considerations can be carried out for sample minima. A given distribution function belongs to the minimum domain of attraction of one of three distributions (Frechet, Weibull or Gumbel). We give the characterization of the minimum domain of attraction of Weibull distribution. A distribution function F belongs to the minimum domain of attraction of Gi.y if and only if i]/ > —00 and 1 -F(rjF + uz) y . ,, ^ n hrn -——-- = z7, for all z > 0 , "^0+ l-F{riF + u) where i]/ inf{z € R : F(z) > 0} is the left endpoint of F. TAIL BEHAVIOUR OF SHAPE AND SIZE PARAMETERS In this section we show the relation between the maximum domains of attraction of the shape and size parameters of profiles and particles. In the proofs we will often use the following lemma. It is a generalization of Lemma 1.2.1 in de Haan (1975) and it can also be found in Kötzer and Molchanov (2006) as Lemma 2.4. Lemma 1. Let /(•,•) and g(-,-) be positive functions such that both fCO fCO / f{s,t)dt and / g{s,t)dt Jo Jo are finite for some 0) e (0,°°] and for s lsag(s,t)ät = C . Now we are ready to prove the stability of MDA for shape factors of spheroids. We start with prolate case. Theorem 1. Let Si,..., E„ be independent and identically distributed prolate spheroids with isotropic orientation. Assume that the orientation is independent of semiaxes lengths and that the E; are not spheres with positive probability. Then the following assertions hold. • IfFs > 0. First we consider the maximum domain of attraction of Frechet distribution. Applying Eq. 3 and Lemma 1 we find that l-FT(uz) lim-—- ► oo \ zVT+U lim —j=. U^oo ^ I uz oo 1—F$(v) The case of MDA of Weibull distribution can be treated in a similar way: I-FT(CO-uz) lim o+ \-Ft{(0-u) a -uz r(o s/l+m- - f Ja 1 -Fs(s) -uz 2s3/2 yjs-co+uz dï lim , ,, u^0+ m-u ra 1 -Fs{s) , v/l+ffl-H 2s3/2 s/s-a+u ,_ f » l-Fs(a-vz) , d (œ-nz)\/l + o}-n J0 (co-vz)3/2^/íü^Yz" o+ (œ-u)VT+œ^î ' f" '-ffi-^Uv J" (Vi) ri- , n r lim A/5 lim ( v + zb{v) is strictly increasing on [mo, Substituting .v = v + zb{v) we get lim + u+zb(u) ro> i-rS(s) , yjl+u+zb(u) Ju+zb(u) 2syiyjs-u-zb(u) 1 -Fs(s) lim U—C0- II . (0 1 -Fs(s) y/l+u (u + zb(u))\J\ + u »-"ffl- 11^/1 +u +zb(u) fg-Um) 1 -Fs(v+zb(v)) d s lim f1 Ju (v+zb(v))3/2 ^/v+zb(v)-u-zb(u) (l + zb'(v))dv rCO 1 -Fs(v) The properties of b (u+zb(u))s/\+u _ ii\J \+u+zb(u) 3/? and Eq. 8 imply IMm^co-il+zb'iu)) = 1, 1 -Fs(ii+zb(n)) _ it^co- 1-Fs{u) Applying the mean value theorem and Lemma 1, we conclude with lim,,^- 1 = e □ It means that the distribution function of profile shapes belongs to the same maximum domain of attraction as the distribution function of particle shapes, only the parameter y may differ. Similar result holds for the population of oblate spheroids. Theorem 2. Let Ei,...,E„ be independent and identically distributed oblate spheroids with isotropic orientation. Assume that the orientation is independent of semiaxes lengths and that the E; are not spheres with positive probability. Then the following assertions hold. • IfFs 0. Applying Eq. 6 we can write Proof. The proof proceeds along the same lines as that of Theorem 1, except that we use Eq. 5. We just show how this works for the maximum domain of attraction of Frechet distribution: 1 -fd\c.{uz\ c) lim-—-—,—- oo 1-Fd\c(u I c) c(lt2z2-C2) r oo UZ Jl "Z (a2_c2)3/2 v/a2_„2.2 da i^i c{u2-c2) roc a(l-FA\c(a\c)) " J" (a2—c2)^i2\fi da lim <2—ll2 r oo vz{l-F^c{vz\ca)) U2Z2-C2 j" (l,2-2_c2)3/2.v/^2 zdv »z(u2 -c2) r~ H1-FaICW) dv J" (V2_c2)3/2.v/v2_„2 1-Fa\c{uz\C) lim-—--:—- 1 —FA\C\U I c) = z-(y+1). □ When considering oblate spheroids, the situation is reversed compared to the prolate case. We observe major semiaxis directly and minor semiaxis of the section profile is larger than minor semiaxis of the spheroid. It means that large minor semiaxis can be observed even if spheroid minor semiaxis is small. On the other hand, small minor semiaxis can only FD\A{uz\a) lim ——--:—- Fd\a(u | a) lim "^f'fp , ,„Fcu(c | a)dc u"0+ r"_ QC(°2-"2)_f , Jr\ a) dr Jo r~>—it 2 2i3/2 c\A\y I ") u v u- — c- (a-—c~)D! ~ J" . FclA(vz I a)zdv J0 -2 v/„2_v2(r2 - -2 -2 "i3 /2 I lim )3/2 0+ r«« FcIa (V | a) dv J0 Vu^v2^2-^)3/2 C|AV ' ; which completes the proof. □ We also consider the maximum domain of attraction of Weibull distribution. Theorem 5. Let Ei,..., E„ be independent and identically distributed oblate spheroids with isotropic orientation. Assume that the orientation is independent of semiaxes lengths A and C. LetFc |A(0 | a) = 0. IfFc |A belongs to the maximum domain of attraction of G2>y with y > 1, then FD\A belongs to the maximum domain of attraction of (h a- Proof. We have to show that 1 —FDiA(a — ud | a) lim o+ 1 — FDiA(a — u | a) d. From Eq. 6, we get 1 —Fd\a(ci — u I a) = a(2a-u) fa-u c(1-Fc|a(c|o)) a-it Jo sj{a-u)2-c2{a2-c2)3!2 d c. The assumption £ MI)A(C/2.y) ensures that 1 — Fc\a(c I a) = ia ~ c)rL(a — c) for some slowly varying function L, i.e., liin„ _o L(uz.)/E(u) = 1 for any z > 0. Hence, we may apply the dominated convergence theorem and deduce that lim H—>-0+ 1 — FDu(a — u I a) u ra c(I —Fc\a(c | «)) = 2aJo (a2 — c2)2 dC<°° Therefore, we can conclude that 1 — FDu(a — ud | a) lim o+ 1 — Fd\a(ci — u I a) 1 - Fd\a(ú - ud I a) u d lim --------:—- ud 1 — FD\A{a — u I a) d. □ WORKED EXAMPLE In this section we demonstrate the utility of our theoretical results on the following example. We consider a population of oblate spheroids and assume that the shape factor S has distribution function h's with finite right endpoint 0 < a> < °° and power law behaviour at co. Specifically, 1 -Fs(s)=K(co-s)c 0 Kw-K-1/* 0. This ensures that Fs belongs to the Weibull maximum domain of attraction MI)A(C/2.y). By Theorem 2, distribution function Ft of profile shape factor T belongs to MI')Aff72 r 1,2). Using Eq. 4 we show that Fj is tail-equivalent to a distribution function with power law behaviour at co. Lemma 2. Let Fs satisfy Eq. 9 and Fj be given by Eq. 4, then lim f—>(q- 1 -FT(t) K f 1 - =-B v+ 1 - (co-t)n 1/2 2V® V 2 where B(-,-) is the beta function. Proof. First we rewrite the complementary distribution function of Fj using substitution v = (co — s)/(co — t), tK{(o-Sy l-Fr(t) _ (m-t)y+1/2 Jt 2s3/2V^t{co-t)r+1/2 f1 tKyr(l-y)-1/2 d s Jo 2{co-{co-t)y)3/2 The statement of the lemma now follows by noting that lim and t^co- (o - (o - t)y)3/2 i/ffl f1 yY(l — v)_1^2dv = B(y+ 1,1/2) . Jo □ Recall that the normalizing constants for Fs can be chosen such that cn = co — Fs (1 — l/n) and dn = co. By Eq. 9, we get cn = (nK)~llv, see also Embrechts et a I. (1997), Example 3.3.16. Lemma 2 yields that the possible choice of normalizing constants for Ft is cn = (riK) and dn = co, where y y I 1/2 and K = ^=B(f+ 1/2,1/2). Since we want to estimate cn from estimate of cn, the following relation turns out to be useful, = [Cn 20S (10) In what follows, we consider a particular example that leads to the shape factor with power law behaviour at the right endpoint. A typical particle So is assumed to be oblate spheroid with semiaxes lengths A and C having joint density function for 1 < c < <7 < 2 , fA,c(a,c) = { C (11) 0, otherwise. After straightforward calculation one obtains distribution function Fs of the shape factor S = A2/C2 - 1, Fs(s) = 1-^(3-s)2, 0 < s < 3 . (12) We see that it is of the form Eq. 9 with y= 2, K = 1/9 and co = 3. It is also possible to express Ft by evaluating integral in Eq. 4, Profile shapes + —77;— log - 18 yft 0 < r < 3 . (13) We generate a sample of n = 250 independent isotropic oblate spheroids with semiaxes lengths distributed according to Eq. 11. For each spheroid a local section through its centre is performed, resulting in the collection of ellipses (profiles). Let S\,...,Sn be the shape factors of simulated particles and let Ti,...,Tnbe the shape factors of their profiles. Fig. 2 shows the histogram of particle shape factors together with theoretical density function. Of course, both are not available when observing only sectional data. On the other hand, Fig. 3 shows the histogram of shape factors obtained from central sections, this is the information that we have in practice when dealing with real data. For comparison, the theoretical density function (in practice unknown) of profile shapes is also depicted. Particle shapes 0.0 0.5 1.0 1.5 2.0 2.5 Fig. 3. Histogram of profile shape factors in simulated population of oblate spheroids together with theoretical density function derived from the distribution function given by Eq. 13. Denote by Mn = max(Si,... ,Sn) the maximal shape parameter. Our aim is to predict the distribution of Mn based on observations of local section shapes 7'i..... Tn. We order them from the largest to the smallest, i.e., T^ > T(2) > > T(n) is the order statistics. We suppose that Fs € MDA(G2,r) and that the right endpoint co = 3 is known. By Theorem 2, Ft s