Buffon's Needle 3D

This section is concerned with Buffon's needle problem extended to three dimensions. Instead of an array of lines (or joints), consider an array of parallel planes, a distance h apart. In ³, these planes can be described by z = nh (n = 0, ±1, ±2, …). A needle, length l (l < h), is thrown into the array of planes. What is the probability, P, that the needle transects one of the planes?

Figure 7: Buffon's needle problem extended to three dimensions.

Let (P, Q, R) be the coordinates of the centre of the needle. Let the orientation of the needle be determined by the two angles, F and Q. Denote the distance from the needle's centre and the plane beneath it as Z. The needle is assumed to be IUR in ³. Therefore Z ΠUR[0, h), F ΠUR[0, p) and Q is sine weighted in [0, p/2). Z has density function f_Z(z) = 1/h. F has density function f_F(f) = 1/p. Q has density function f_Q(q) = sinq. Z, F and Q are independent variables and so f_{z, F, Q}(z, f, q) = f_Z(z) f_F(f) f_Q(q). Thus, the triple Z, F, Q has joint density function f(z, f, q) = (sinq)/ph.

Figure 7  shows the needle transects a plane if and only if z £ (l/2)cosq or z ³ h - (l/2)cosq. Clearly, a transect occurring is not dependent upon f. Let B₃ = {(z, f, q) : z £ (l/2)cosq or z ³ h - (l/2)cosq}. A transect occurs for all (Z, F, Q) Î B₃. A trivial extension of RandomVar 1-(11) leads to the following:

    (14)

The following work assumes h is known with l to be estimated. The needle transects a plane with probability l/2h and falls between planes with probability 1 - l/2h. Consider the discrete random variable Q that maps the outcomes needle transects a plane and needle falls between planes to the real numbers 1 and 0, respectively. According to RandomVar 1-(1), Q has expected value EQ = l/2h and therefore

   (15)

An estimate of l can be obtained from a single throw of the needle:

   (16)

Curve Length 3D

Buffon's needle in three dimensions can be extended to estimate the length, L, of a curve, C, in ³. The curve C is thrown isotropic uniform randomly into the array of planes.

Suppose the curve C is approximated by a union of line segments:

If the length of y_i is l_i (1 £ i £ n), an approximation of the curve length of C is:

By comparison with (14), the probability that the line segment y_i intersects a joint is l_i/2h. Furthermore, (15) implies the length of y_i is l_i = 2 h EQ_i. Q_i is a discrete random variable that maps the outcomes y_i transects a plane and y_i falls between planes to the real numbers 1 and 0, respectively. Let the total number of transects, Q, between the planes and Y be denoted as the sum Q₁ + Q₂ + … + Q_n. Taking into account all line segments:

   (17)

An estimate of L can be obtained from a single throw of Y. An estimate of L is

   (18)

Finally, consider again the systematic set of parallel planes a distance h apart. The average area of test plane likely to fall within a certain volume can be calculated. Figure 8(a) shows an area of test plane, A, and the volume it occupies. The test system of planes has area per unit volume A/V = A/Ah = 1/h. Next, substitute h = V/A into (17). The result is an important formula originally due to Saltykov (1946):

   (19)

An estimate of L can be obtained from a single throw of Y:

   (20)

Suppose a needle, length L (L > h), is thrown with IUR position onto a test system of parallel planes. is greatest when the needle is perpendicular to the test planes and smallest when the needle is parallel to the test planes. A better estimation of L is obtained if the test system of parallel planes is replaced by a test system of three sets of mutually orthogonal planes (see Figure 8(b)). Such a test system has three times as much area per unit volume as the parallel plane test system so that A/V = 3×(1/h).

(a)                                               (b)

Figure 8: Test systems of (a) parallel planes (A/V = 1/h) and (b) mutually orthogonal planes (A/V = 3/h). The figures shows each test system with the volume, V, an area of test surface, A, occupies.

Surface Density

Equation (19) says that the length of C per unit volume is equal to twice the number of transects between C and the test system, per unit area of test system. The dual situation can be conceived if the test system is replaced by a stationary, bounded object in ³, M, whose surface area is S and whose volume is V. Suppose S and V are unknown. L/V is then the length of C lying within M, while EQ/A becomes the number of intersections between C and M, per unit area of M. By an additional change of notation of A into S and of Q (transect of C with a plane) into I (intersection produced by C in a surface) the following is obtained:

   (21)

Equation (21) determining the surface density, S/V, of M is another fundamental equation of stereology originally obtained by Saltykov (1945). In practise, the (test) curve C is replaced by an IUR test system of test lines. I and L are random variables. I denotes the number of intersections between the object surface and the test lines, L the length of test lines falling within the object.

An IUR line in ³ that passes through M can be generated as follows. Let u₁ Î UR[0, 1), u₂ Î UR[0, 1), u₃ Î UR[0, 1) and u₄ Î UR[0, 1). Construct a line or isotropic axis with orientation q = cos^-1(1 - 2u₁) and f = pu₂. As shown in Figure 9, fix a plane perpendicular to the isotropic axis and project the object, M, onto that plane. Draw a square on the plane (side length l), which contains the projection of M. From a point that has UR position within the square, namely (u₃×l, u₄×l), draw a second axis perpendicular to the plane. If this axis hits M (i.e. the point (u₃×l, u₄×l) is contained within the projection of M), then it can be considered an IUR line in ³ that passes through M. Otherwise the procedure must be restarted with 4 new random numbers.

Figure 9: The construction of an IUR line.

Suppose a UR systematic set of points, a distance u apart, are positioned along each IUR test line. The length of test line associated with each point, l/p, is u. If the expected number of points falling within the object is EP, the length of L is EP×(l/p) (see Figure 10). Therefore, an estimate of S/V can be obtained from a single throw of the IUR test line system:

   (22)

Figure 10: A method for estimating the length, L, of test line falling within an object. In this example, the object is a sphere and an unbiased estimate of L is 3u.

Surface Density - IUR Sections

In practice, surface density estimation is achieved by taking IUR sections through the object of interest, M, and superimposing an IUR test system of lines onto each section. Subsequently, surface area may be calculated by multiplying the estimate of surface density by the reference volume, which may be obtained using the Cavalieri method.

An IUR section (or plane) in ³ that cuts M is generated by the following procedure. Let u₁ Î UR[0, 1), u₂ Î UR[0, 1) and u₃ Î UR[0, 1) so that q = cos^-1(1 - 2u₁) and f = 2pu₂ define an isotropic orientation, n Π³. Let h be the radius of a sphere that encompasses M. As shown in Figure 11 the point a Π³, found by moving a distance h×u₃ from the sphere's centre in the direction n, and the orientation n together define an IUR plane in ³ that may cut M. If the plane misses M, the procedure must be restarted with 3 new random numbers.

Figure 11: The construction of an IUR section.

When working with biological material it can be extremely cumbersome to generate IUR, physical sections. Mattfeldt et al. describe a way to generate IUR sections using the "orientator" (1990). Nyengaard and Gundersen (1992) describe an alternative approach using the "isector". However, non-invasive scanning techniques such as MRI enable a complete rendering of 3D structures, from which image analysis packages can generate IUR sections.

Given a set of IUR sections through M, it is relatively easy to superimpose IUR test systems of lines onto each section. Generally, it is more efficient to superimpose IUR test systems comprising parallel or orthogonal lines. An IUR test system of orthogonal lines, a distance h apart, can be generated as follows. Let u₁ Î UR[0, 1), u₂ Î UR[0, 1), u₃ Î UR[0, 1) and q = 2pu₁. As shown in Figure 12(a), O Π² is a fixed origin on the IUR section and q is an angle relative to the fixed axis, O_x. The points a_i, b_i Π² (-¥ < i < ¥) are found by moving, respectively, distances h×u₂ + i×h and h×u₃ + i×h from O in the directions q and q + p/2. An IUR test system of orthogonal lines is constructed by drawing two sets of lines with orientations q + p/2 and q through the points a_i and b_i.

Figure 12(b) shows an IUR section through an object M with an IUR test system of orthogonal lines overlain. The orthogonal grid of lines has 5 of its vertices within M. The length of line, l, and area, a, associated with each vertex of the grid is l = 2h and a = h². Therefore, an estimate of the length, L, of test line falling within M is 5×2h = 10h. The number of intersections, I = 10. The test system of orthogonal lines is said to have length per unit area l/a = 2/h. A similar test system comprising a single set of parallel lines has length per unit area l/a = h/h² = 1/h.

(a)                                                   (b)

Figure 12: (a) The construction of an IUR test system of orthogonal lines. (b) The calculation of length per unit area, l/a, for a test system of orthogonal lines.

Test systems of parallel and orthogonal lines are made up of test line elements that may be rearranged into any arbitrary curve shape whose mean length per unit area, l/a, is known. Therefore, an alternative method is to overlay UR test systems of rolling circles (similar to that shown in Figure 6) on IUR sections through M and count intersections between the profile of M and the rolling circles. An estimate of the length, L, of test line falling within M is obtained by counting points falling within M that are placed at quarter circle (pr/2) intervals along each rolling circle.

Surface Area - Vertical Sections

It is often more convenient to estimate surface area from vertical sections. A vertical section, as shown in Figure 13, is a section that is perpendicular to an arbitrary chosen horizontal plane and its intersection with the horizontal plane is an IUR line on the horizontal plane. On anatomical datasets, coronal, sagittal and sagittal-oblique sections can be obtained as vertical sections when the horizontal plane is an axial section.

Figure 13: Vertical sections.

Any straight line in ³ can be contained in a vertical section. Therefore, any stereological procedure requiring IUR test lines, such as (22), may be implemented using vertical sections as supports for the lines. A test line on a vertical section is effectively IUR in ³ if the angle q between the test line and the vertical direction has probability density function sinq (see Random Position and Orientation) and, given q, the position of the test line is uniform random on any bounded interval perpendicular to the test line. Thus, a uniform random sine-weighted line on a vertical section is effectively IUR in ³. Instead of repeated sampling of sine-weighted straight lines on a vertical section it is convenient to use a test curve that encapsulates the sine weighting property for all q (0 £ q < p), in the sense that the length of element of arc ds(q) whose tangent makes an angle q with the vertical, is proportional to sinq. The cycloid (illustrated in Figure 14(a)), aligned with its minor axis parallel to the vertical direction, of classical geometry is such a curve (Baddeley et al., 1986).

(a)                                                  (b)

Figure 14: (a) A cycloid, parametrically defined as: x(q) = (q - sinq)r and y(q) = (1 – cosq)r with 0 £ q £ p. (b) A test system of rolling cycloids or cycloid chains (Cruz-Orive and Howard, 1995).

Given vertical sections through a structure with UR cycloids overlain, (22) can be used to estimate the surface density. I denotes the number of intersections between the boundary of the surface of interest on the sections and the cycloids and L the length of cycloid within the structure. Suppose a point is marked at the beginning of each cycloid. The length of cycloid associated with each point, l/p, is 4r. If the expected number of points falling within the object is EP, the length of L is EP×(l/p).

A better estimate of the surface density is achieved if the cycloids are systematically placed on the vertical sections. Figure 14(b) shows a UR test system of rolling cycloids a distance 4r apart. A rolling cycloid is defined as

   (23)

   (24)

where [x] means the integer part of x. The test system can be thought of as a tessellation constructed from a mosaic of fundamental tiles, J₀. J₀ occupies an area a = 4pr.4r and contains a length of rolling cycloid, l = 16r. Therefore, the test system has length per unit area, l/a = 16r/(4pr.4r) = 1/pr.

Surface Area - Exhuastive Vertical Sections

In the case where several series of vertical sections, separated by a distance T, sample the object of interest M from one end to the other (i.e. exhaustively), reference volume can be determined from the same images as provide the estimate of surface density. If the Cavalieri estimate of volume (5) is substituted into (21), an equation determining surface area directly by intersection counting is obtained:

   (25)

The best estimates of surface area are achieved if M is sampled systematically. In Figure 15, M is sampled systematically by n = 4 series of vertical sections. Series i (i = 0...n - 1) has orientation fi = f₀ + i×(p/n) where f₀ Î UR[0, p/n). The orientations are shown as dotted lines emanating from an arbitrarily chosen point P on the horizontal plane. In Figure 15, P is positioned near the centre of the projection of M onto the horizontal plane. The offsets, z_i, for each of the n series of vertical sections are chosen at random without replacement from the set of values (j/n)×T (0 £ j < n).

Figure 15: Illustration of a method for obtaining exhaustive vertical sections through a bounded object, M. The view is from above M, looking down on the horizontal plane. In this example, the number of series of vertical sections, n, is 4.

Each series of vertical sections through M, with UR test systems of rolling cycloids overlain, contributes an estimate of surface area, S_i, according to (25). A final estimate of object surface area, S, is calculated from the mean of the surface area estimates obtained from each series so that S = (S₀ + … + S_n)/_n. This result was first recognised by Michel and Cruz-Orive (1988).

References

BADDELEY, A. J., GUNDERSEN, H. J. G. and CRUZ-ORIVE, L. M. Estimation of surface area from vertical sections. J. Microsc., 142, 259-276 (1986).

CRUZ-ORIVE, L. M. and HOWARD, C. V. Estimation of individual feature surface area with the vertical spatial grid. J. Microsc., 178, 146-151 (1995).

MATTFELDT, T., MALL, G., GHAREHBAGHI, H. and MÖLLER, P. Estimation of surface area and length with the orientator. J. Microsc., 159, 301-317 (1990).

MICHEL, R. P. and CRUZ-ORIVE, L. M. Application of the Cavalieri principle and vertical sections method to lung: estimation of volume and pleural surface area. J. Microsc., 150, 117-136 (1988).

NYENGAARD, J. R. and GUNDERSEN, H. J. G. The isector: a simple and direct method for generating isotropic, uniform random sections from small specimens. J. Microsc., 165, 427-431 (1992).

SALTYKOV, S. A. Stereometric Metallography, 1st edn. (In Russian). State Publishing House for Metals Sciences, Moscow (1945).

SALTYKOV, S. A. The method of intersections in metallography (In Russian). Zavodskaja laboratorija, 12, 816-825 (1946).