Sedimentological analysis of grain fabric has paid scant attention to grain shape. However, the information of grain orientation is inseparable from that of shape in three-dimensional fabric analysis. Not only should the dominant major-axis orientations be recognized, but so should the dominant combinations of shapes and orientations of grains. The aim of this paper is to demonstrate that such combinations can be identified by density-based cluster analysis in a five-dimensional parameter space, where a point represents a specific combination of the shape and orientation of a grain approximated by a triaxial ellipsoid. We tested the present method using an artificial data set. The data were successfully classified into correct groups. Next, we applied it to a data set obtained by X-ray computed microtomography from some 5000 sand grains deposited in an experimental flume. We show that triaxial grains, which have principal axes with distinctive radii, have major axes in the paleocurrent orientation and roller-shaped grains in the transverse orientation. Both grain types have vertical minor axes. Those orientations are not preconditioned for statistical analysis but are recognized as significant; this suggests the potential of the method in three-dimensional fabric analysis for applications in sedimentology. The method can be applied to the statistical processing of ellipsoidal objects including, for example, deformed grains, stress ellipsoids, and magnetic susceptibility ellipsoids.
The three-dimensional fabric of rocks is a clue to the solutions of various problems including those of tectonic deformations (Ketcham, 2005b; Shao and Wang, 1984) and environments. For example, the fabric affects anisotropic fluid flows and associated chemical processes at depths (e.g., Phillips, 1990). A memorable example occurred when paleocurrent direction was inferred from the stereo images of imbricated rocks near the Mars Pathfinder landing site (e.g., Golombek et al., 1999; Basilevsky et al., 1999).
The recent development of X-ray computed microtomography has enabled us to capture the three-dimensional shapes and orientations of grains in rocks (e.g., Mees et al., 2003). However, most studies on grain fabric use the traditional methodology developed when grains were measured on two-dimensional surfaces, e.g., thin sections, outcrops, etc. Specifically, major-axis orientations are often used in those studies.
Here, we present a statistical technique that can be used to recognize dominant combinations of orientations and shapes of grains from three-dimensional data, where the shapes are approximated by ellipsoids. The technique is applicable not only to sedimentary grains but also to the statistical processing of ellipsoidal objects, including deformed grains in metamorphic rocks (e.g., Ikeda et al., 2000; Ketcham, 2005b) and stress and magnetic susceptibility ellipsoids (e.g., Orife and Lisle, 2003; Jelinek, 1981).
The term “fabric” in sedimentology refers to the packing and orientations of grains (Martini, 1978). This study is concerned not with packing but only with orientation. Sedimentary grains characteristically have preferred orientations to some degree (Griffiths, 1967), and the orientations have been often related to paleocurrent directions (Allen, 1982, chapter 5). Most types of deposits show elongate grains oriented in the flow direction (e.g., Middleton, 2003). Recent studies have demonstrated that paleohydraulic conditions can be distinguished using grain fabric even if the conditions resulted in indistinguishable sedimentary structures (Endo et al., 2002; Yagishita et al., 2004). The analysis of sedimentary grain fabric has been, accordingly, concerned with the global maximum or local maxima in the frequency distribution of major-axis orientations of grains. The shapes of grains or clasts have been utilized to assess the paleoenvironment under which they were deposited (e.g., Cailleux, 1945; Dobkins and Folk, 1970; Gale, 1990; Howard, 1992; Benn, 1994a; Benn and Ringrose, 2001; Oakey et al., 2005), though these studies did not include orientation data.
The difficulty of three-dimensional fabric analysis comes primarily from the complexity of three-dimensional data due to the inseparable shapes and orientations of grains, even if grains are simplistically approximated by ellipsoids. Orientations of elongated or platy grains are more important than those of nearly spherical ones for the study of the physics of sedimentation (Middleton and Southard, 1977; Allen, 1982, p. 179). However, not only aspect ratios but also shape factors require attention (Sneed and Folk, 1958): the shape factor describes the spectrum between prolate and oblate ellipsoids. Namely, major axes are more important than intermediate and minor axes for elongated grains (Fig. 1), but the minor axes are significant for oblate grains because they designate grain attitudes. A triaxial ellipsoid, which has three different principal radii, has major- and minor-axis orientations that are equally important in sedimentation (Zapryanov and Tabakova, 1999). Therefore, conventional fabric analysis in sedimentology that investigates only major-axis orientations or neglects disc-like grains is ineffectual for three-dimensional data. The orientations of grains should be studied jointly with the shapes of grains.
Recent advancement in X-ray computed tomography has enabled us to recognize three-dimensional shapes of fine grains in a rock (e.g., Mees et al., 2003; Ketcham, 2005a). Such technology provides a massive amount of data from thousands of grains, giving rise to another difficulty—the way to summarize such a massive and complicated data set and to extract significant information. Three-dimensional fabric analysis is expected to recognize statistically dominant and possibly mechanically significant combinations of shapes and orientations of grains from a large data set.
We meet this demand by adapting cluster analysis to search for extremal points in the joint frequency distribution of shapes and orientations of grains. In this article, we begin by introducing a five-dimensional parameter space: the data are transformed to points in the space, to which clustering is applied in order to recognize significant extrema of the joint distribution. Second, we introduce a clustering technique. Third, an artificial data set is used to test our technique. And, finally, data obtained through X-ray computed tomography from a sand layer are used to show that the present method has potential for sedimentological studies. The present method is an extension of that of Yamaji and Masuda (2005), who utilized an Elliott plot for two-dimensional fabric analysis (Fig. 2). This is the first attempt to apply a modern statistical technique to the understanding of three-dimensional grain fabric. It is not the purpose of this paper to investigate the way grains became oriented.
The shapes and orientations of grains are parameterized as follows. Throughout this study, grains are approximated by ellipsoids, and their sizes are neglected (Sneed and Folk, 1958). We always deal with the shape and orientation of a grain as grouped data.
For the case of two-dimensional fabric analysis, we deal with the aspect ratio R and major axis orientation φ of an elliptical grain. The grouped data of an ellipse have a one-to-one correspondence to a point on a plane with the two-dimensional Cartesian coordinates O-xpxq, where xp = (logeR)cos 2φ and xq = (logeR)sin 2φ (Fig. 2). The origin (O) of this two-dimensional parameter space represents circular grains, whereas points distant from the origin stand for elongated grains. The two points (xp, xq) and (−xp, −xq) represent grains that have the same aspect ratios and major-axis orientations perpendicular to each other. Suppose that several groups of grains are given, where grains are classified by their shapes and orientations. The points corresponding to the grains of a group make a cluster on the plane. Also, the cluster center stands for the representative shape and orientation for the members of the group. Accordingly, recognition of the joint frequency of shape and orientation results in the cluster analysis of the points on the plane (Elliott, 1970; Yamaji and Masuda, 2005). The plane works as a convenient parameter space for two-dimensional grain fabric.
Unfortunately, the parameter space for three-dimensional grain fabric is inevitably high-dimensional and difficult to visualize. We explain the necessary conditions for the space here: the full definition of the space is given in Appendix A.
Unfortunately, the simple parameter space that is spanned by the axes of the Cartesian coordinates O-RΦsdp is not appropriate for our purpose, where O stands for the origin of the coordinates, because a point in the space does not have a one-to-one correspondence to the grouped data. First, a triaxial ellipsoid has orthorhombic symmetry, resulting in the periodicity of 180° rotation about any rotation axis. Accordingly, countless points in the space represent the same figure in the physical space. Second, p is not significant for describing the attitudes of nearly prolate grains (Fig. 1). For those grains, p can have different values even if they have analogous shapes and attitudes. To make matters worse, the angles, p, d, and s, cannot be defined for oblate grains because they have indistinguishable major and intermediate axes.
Accordingly, we adapt the elaborate five-dimensional parameter space of Sato and Yamaji (2006), who modified the parameter space proposed in plasticity theory by Ilyushin (1963). Sato and Yamaji used the space for statistical processing of stress ellipsoids. Design of a parameter space is critical to our problem, but it is somewhat complicated. So, the details of the parameter space are given in Appendix A. The space is designed to overcome the problem of axially symmetric grains whose principal orientations, except for those coinciding with the symmetry axes, are uncertain. Namely, a point in the space has a one-to-one correspondence to the combination of the shape and orientation of an ellipsoidal grain, even if the grain has an axial symmetry.
Let O-x1x2x3x4x5 be the Cartesian coordinates of the five-dimensional parameter space, and x = (x1, x2, x3, x4, x5)T be a position vector. O is the origin of this coordinate system and represents spherical grains. In what follows, a position vector is used to indicate the point at the end of the vector, and vice versa. Given N grains, fabric data are transformed to N data points, x(1),…, x(N), in the space. Cluster analysis is applied to those points.
Dozens of clustering techniques have been proposed (Tan et al., 2006). We employ density-based multiscale clustering that searches for high-density regions of data points in the parameter space because this method does not need the number of clusters a priori and is tolerant to background distribution of points between clusters (Hinneburg and Keim, 1998; Nakamura and Kehtarnavaz, 1998).
The continuous function D(x, i) has its maximum at the position x(i) and decreases monotonously with increasing distance from x(i). The scale size defines the resolution of a point (Fig. 3). If σ is smaller than the spacing between points, each point is resolved as a spike in the function D(x). On the other hand, very large σ results in a single maximum that represents the entire point distribution.
If there are clusters in the parameter space, they are identified as peaks of D(x) with appropriate values of σ, but increasing σ amalgamates neighboring peaks to form larger ones (Fig. 3). Consequently, the number of peaks changes with increasing value of this parameter, but the peaks corresponding to statistically significant clusters persist for a wide range of σ (Hinneburg and Keim, 1998; Nakamura and Kehtarnavaz, 1998). This persistence is the basis for identifying statistically significant peaks.
The density should be corrected to account for the geometry of the parameter space. It is obvious in Figure 4 that the density is inevitably greater in the vicinity of the origin than in the areas distant from the origin, even if nearly spherical grains are rare. The high density near the origin is an artifact that results from the assignment of an aspect ratio to the distance from the origin (Appendix A). This must be corrected.
Detection of Clusters
Peaks of the function F(x) were detected by a hill-climbing technique starting from data points. Specifically, we used the simplex method (Nocedal and Wright, 2000) for this purpose. Once a peak was found at x, the pair of shape and orientation represented by the peak was determined (Appendix A). Climbers from different data points in a cluster eventually arrive at the same peak. So, the ellipsoid corresponding to the position of the peak is thought to represent the ellipsoids corresponding to those data points.
The present method was tested using an artificial data set with a known answer. The data included 150 ellipsoidal grains and were subdivided into three groups A, B, and C (Fig. 6A). Group A included 70 nearly prolate grains, the major and minor axes of which were scattered around the vertical and north-south, respectively. The grains had mean Φ and R at 0.10 and 2.00, respectively (Figs. 6B and 6C). Approximately oblate grains made up Group B. Their mean Φ and R were 0.90 and 1.67, respectively. The major and minor axes were scattered around the vertical and an east-west trend, respectively. Group C was composed of triaxial grains with the mean Φ and R at 0.50 and 2.5, respectively. The mean major-axis orientation was inclined at 60°, whereas the mean minor-axis orientation was inclined at 28° (Fig. 6A).
Groups A and B were, therefore, indistinguishable from their major-axis orientations alone. Their minor-axis orientations are separated by ∼90°, but the minor axes of Group C filled the gap between them. As a result, the three groups were difficult to distinguish from their principal orientations. Although the frequency distribution of R had three modes, which corresponded to the groups, the distributions of the groups had overlapping tails (Fig. 6B). On the other hand, the distributions of Φ were clearly separated (Fig. 6C).
The points corresponding to the data make three distinctive clusters in five-dimensional parameter space. Figure 7 shows the three-dimensional plots of the points (x1(i), x2(i), x3(i)) and (x2(i), x4(i), x5(i)), where xj(i) is the jth component of the ith member of x(1),…, x(N), and the data points correspond to the 150 grains. The plots show three distinctive clusters.
Figure 8 shows the results of the method with various scale sizes: the groups were successfully identified from the data. Groups A, B, and C were recognized as the long-standing major clusters for increasing σ. The clusters corresponding to Groups A and B were distinguished for smaller σ values, but the clusters were amalgamated for larger σ values due to the small Euclidean distance between the cluster centers (Fig. 7C). Group C was detected with larger σ values because the cluster standing for the group was more scattered than the other clusters (Fig. 7). Figure 6D depicts the three dominant combinations of shapes and attitudes of the grains.
It may seem nonintuitive that the peaks corresponding to Groups A and B are combined but are distinctive from the peak corresponding to Group C when σ is high (Fig. 8), even though Group B has a very well-defined minor axis that is 90° away from the minor axis of Group A (Fig. 6A). However, the amalgamation of the clusters of Groups A and B is a natural consequence, since the groups have extreme values, i.e., Φ ≈ 0 and 1, respectively. On one hand, grains with Φ ≈ 0 have well-defined major-axis orientations, but their minor-axis orientations are scattered on the plane perpendicular to the cluster center of the major-axis orientations. On the other hand, the opposite is true for grains with Φ ≈ 1. Therefore, the 90° difference between the cluster centers of the minor-axis orientations of Groups A and B is not as significant as the differences between Groups A and C and between B and C. The distances between the cluster centers in our five-dimensional parameter space corresponding to the three groups are shown in Figure 7C. The distance between Groups A and B is 0.602, but Group C is separated by distances of 0.787 and 0.931 from Groups A and B, respectively.
APPLICATION TO SAND IN AN EXPERIMANTAL FLUME
The present method was applied to a data set obtained from a planar bed of sand deposited under a unidirectional flow in an experimental flume. The layer was solidified with resin, and an oriented sample was taken from the layer. Grains in the sample were analyzed with X-ray computed microtomography at Japan Synchrotron Radiation Research Institute. The shapes and orientations of 4549 grains were determined through the procedure given in Yokokawa et al. (2004). The spatial resolution of the tomographic data was ∼6 μm, and the mean grain size was 100 ±40 μm, where the geometric mean of principal diameters, 2(a1a2a3)1/3, was used as the size of a grain. The grains had a positive skewness of ∼0.32 in the frequency distribution of Φ values: nearly oblate grains were slightly fewer than nearly prolate ones (Fig. 9). R and logeR of the grains had modes at ∼2 and ∼0.6, respectively.
The stereoplots in Figure 10 show the data. Because of the large number of grains, we do not use a “tadpole plot” as in Figure 6A. The grains with R values smaller than 1.5 appear to have random orientations. In contrast, there are bluish and greenish nebulas in the middle and lower panels of Figure 10, indicating that prolate and triaxial grains with R ≥ 1.5 had preferred orientations. It is obvious from the figure that there are many grains with nearly vertical minor axes. Major axes with large inclinations are fewer than those with low to intermediate inclinations.
However, it should be stressed that such concentrations of principal orientations do not necessarily indicate extrema of the joint distribution of shape and orientations of grains but merely show the concentration of either major- or minor-axis orientations. Grains that belong to a cluster of the former axes could have different minor-axis orientations or different aspect ratios. It is difficult to recognize the joint distribution using illustrations due to the high dimensionality of the parameter space.
The result of our method applied to the data is shown in Figure 11. Two groups of grains were recognized as long-standing local maxima of F(x) for seven σ values.
Grains of the most significant group had approximately vertical minor axes and major axes shallowly plunging to the west (upstream), and had representative R and Φ values in the ranges 2.1–2.8 and 0.3–0.4, respectively. Taking their median values at R = 2.45 and Φ = 0.35, we have an ellipsoid with principal radii of a3 = 1.0, a2 ≈ 1.5, and a1 = 2.45 as an example of the triaxial grains. This falls in Sneed and Folk's (1958) “bladed” class of grains.
The second prominent cluster in the parameter space represented nearly prolate grains (Φ = 0.2–0.3), with the major axes perpendicular to the flume current. The grains had representative aspect ratios of 1.9–2.3. Therefore, grains of the second group were represented by an ellipsoid with the principal radii of a3 = 1.0, a2 ≈ 1.25, and a1 = 2.0, which belongs to the class of “elongated” grains of Sneed and Folk (1958).
Based on rough analysis of the sources of error in our method, we found that the uncertainty in R and Φ arose from the migration of the centers of clusters in parameter space, and that of principal axes came from the resolution of the positions of the centers for the detected bladed and elongate classes.
The principal radii and orientations of grains had uncertainty of ∼10%, because the spatial resolution of the tomography (∼6 μm) was about one-tenth of the mean radius of the grains (50 μm). The components of E had similar uncertainty, which propagated to the precision of x roughly with the same magnitude, since the transformation between E and x was linear. The extremal points were detected with the range of σ between 0.10 and 0.22. This lower bound is similar to the precision of x. Therefore, considering the precision and the range of σ, we take the value of 0.2 = 20% as the upper bound of error in the position of peaks in parameter space. Consequently, the detected aspect ratios and shape factors had errors of ∼20%, because the aspect ratios (R = 1.9–2.8) and shape factors (Φ = 0.2–0.4) indicate that the E matrices corresponding to the detected maxima were not ill-conditioned. The principal orientations had errors at ∼0.2 radians ≈ 11°.
These uncertainties are smaller than the ranges of R and Φ that were presented in the previous subsection. For example, the most dominant combination had a range of R between 2.1 and 2.8. These ranges arose from the migration of the extrema points in parameter space by increasing σ, instead of the uncertainty of positions of the extrema themselves, which were measured by σ. Consequently, the ranges indicate the uncertainty of the result. In contrast, the principal orientations of the detected bladed and elongated classes were stable during the increasing σ, so they had an uncertainty of ∼11°.
Major- and minor-axis orientations usually have different confidence regions. Not only do various sources of error affect the regions, but so does the Φ value that corresponds to the detected cluster center. The minor-axis orientation has a large uncertainty compared to the major-axis orientation, and the probable minor-axis orientations make a girdle with the pole parallel to the cluster of the major-axis orientations (Sato and Yamaji, 2006; Fig. 11). The opposite is true for the cases of Φ ≈ 1. Both the orientations have similar confidence regions if Φ = ∼0.5. It is difficult to discuss the confidence region of the Φ value independently from that of principal orientations. The confidence regions of principal orientations and Φ values are addressed by Sato and Yamaji (2006) and Yamaji and Sato (2006) using theoretical analyses and numerical experiments. The confidence regions are represented by the spread of clusters of data points on the unit hypersphere. Namely, the grains with the same principal orientations but with the extreme values of Φ = 0 and 1 have corresponding points on the hypersphere that are separated by an angle of 60° (Yamaji and Sato, 2006, p. 936). If two grains have the same value of Φ = 0, the 35° rotation of the grains corresponds about an angular distance of ∼60° on the hypersphere (Figure 1 of Yamaji and Sato, 2006). If both the grains have a value of Φ = 0.5, the same angular distance corresponds to a 90° rotation of the major axis about the minor axis. However, the theoretical expression of the confidence region of R is not well understood.
Sedimentological studies of grain fabric have been concerned only with major-axis orientations except for a limited number of studies (e.g., Lindsay, 1968). However, information about shapes and orientations is inseparable, and the joint frequency of shape and orientation of grains is important. The present method enables us to recognize their dominant combinations from a great number of grains.
The dominant groups of grains determined here are concordant with the bimodal major-axis orientations that have been often observed in sediments, i.e., the paleocurrent and its transverse orientations. The grains oriented in the transverse orientation were generally roller-shaped. Grains of both groups had nearly vertical minor axes. The three orientations, i.e., the vertical, the paleocurrent, and the transverse orientation, had not been incorporated as special ones in our method beforehand, but appeared as dominant orientations in the result, suggesting that the method has potential for studying three-dimensional grain fabric of sediments.
The present method successfully discriminated the three groups from test data. It should be noted that grouping was not possible either in the major-axis orientations or in the minor axes (Fig. 6A). Both the orientations are needed for the classification. Namely, the joint probability of the major- and minor-axis orientations is important. However, grouping was possible without the present method for this example only by using the three modes of the distribution of Φ values (Fig. 6C). In addition, three clusters are clearly separated in Figure 7. The distribution denotes the projection of data points in five-dimensional space onto one-dimensional space, and the points in Figure 7 are the shadows of the data points. Grouping is possible in the lower-dimensional parameter space if the shadows make clearly separated clusters in the space. However, there are data points that have distinctive clusters in the five-dimensional space but have a continuous distribution in lower-dimensional parameter spaces. Our method can identify such groups.
It is easy to take grain sizes into account, although we do not deal with them in this work. To do this, size is taken as the sixth coordinate, and clustering is conducted in six-dimensional parameter space. However, it should be noted that the definition of size is not unique, depending on the expectation of researchers about the way grains are oriented, which is beyond the scope of this study. For example, if a thin platy grain (a1, a2 ≫ a3), such as a flake of mica, is as important as a larger fat grain (a1 ≈ a2 ≈ a3) in the dynamics of sedimentation, the arithmetic mean of principal diameters, (2a1 + 2a2 + 2a3)/3, is better than their geometric mean, 2(a1a2a3)1/3, to indicate their sizes. The reason for this is that the thin and fat grains have similar arithmetic means, but the thin grain has a much smaller geometric mean than the fat one. In contrast, if one assumes that the volume or mass of a grain is significant to understanding of the formation of grain fabric, the geometric mean is the straightforward measure of grain size, since the volume is proportional to the cube of the geometric mean.
It is one of the advantages of our formulation that the five-dimensional parameter space has no singular point. In addition, every point in the space has one-to-one correspondence to the combination of shape and orientation of an ellipsoid. These advantages are comparable with those of an Elliott plot (Fig. 2) for two-dimensional fabric over the parameter space with the coordinates (logeR, 𝛉), which are often used in structural geology (Ramsay and Huber, 1983). The parameter space with the coordinates (logeR, 𝛉) has singularity: points in the domain R < 1 are meaningless, and entire points along the R-coordinate axis stand for a circle. The one-to-one correspondence does not hold. Regarding three-dimensional fabric, Φ has the closed interval [0, 1]. Therefore, if Φ is taken as one of the coordinates of parameter space, the space has limits at Φ = 0 and 1, and points beyond them are meaningless. Obviously, the hypersphere has no such limit (Fig. A1).
APPENDIX A. PARAMETER SPACE
We construct the parameter space corresponding firstly for the pair of Φ and the orientation of an ellipsoid to have a one-to-one correspondence to a point in the space. The correspondence is briefly explained by Yamaji and Sato (2006). R is embedded in the space afterward.
Shape Factor and Principal Orientations
The Φ value and orientation of an ellipsoidal grain are represented by a point on the unit hypersphere in five-dimensional Euclidean space, which has been used for the statistical processing of stress ellipsoids by Sato and Yamaji (2006).
Thanks to the factor 1/Λ in 201, the vector e has a unit length |e| = 1 (Sato and Yamaji, 2006). If we think of e as a position vector, its end point is plotted upon the surface of the unit hypersphere in a five-dimensional space. As a result, a point on the hypersphere has a one-to-one correspondence to the combination of Φ and principal axes of an ellipsoid. The points that are close to each other on the hyper-sphere represent ellipsoids with analogous pairs of attitudes and shape factors in physical space. An antipodal pair of points stands for ellipsoids with “opposite” shapes and orientations (Yamaji and Sato, 2006): they have interchanged major and minor axes and shape factors of Φ and 1 – Φ (Fig. A1). This correspondence is analogous to that in the two-dimensional case. Namely, the points with positions that are symmetric with respect to the origin of an Elliott plot represent the congruous ellipses with principal axes perpendicular to each other (Fig. 2).
It is known that ellipsoids with principal orientations in common have corresponding points along a great circle on the hypersphere (Yamaji and Sato, 2006). Namely, the e vectors of those ellipsoids are coplanar. It is an important attribute of this hypersphere that although Φ values are limited in the interval [0, 1], a point representing a prolate or oblate ellipsoid is not plotted at the end of our parameter space. The point is on the hypersphere, and it is surrounded by the points that represent ellipsoids with similar shapes or orientations (Fig. A1).
Given a point x in the space, the shape and orientation of the grain that are represented by the point are calculated as follows. The aspect ratio is readily calculated as R = e|x|. Substituting the components of the unit vector e = x/|x| into 05, we obtain the components of E. The principal axes and shape factor of the grain are determined by solving the eigen-problem of this matrix. The eigenvectors of E indicate the principal orientations, and the Φ value is obtained through 03. Finally, the normalized principal radii are calculated via 101.
The authors thank Richard Ketcham and Chris Paola for their constructive reviews and comments that have helped to improve this paper. We acknowledge F. Masuda for the use of his experimental flume and T. Nakano for many pieces of advice on his software for obtaining three-dimensional fabric data through X-ray computed tomography. Thanks are also due to H. Naruse for the collaborative acquisition of the natural data. The synchrotron radiation experiments were performed at the BL20B2 in the SPring-8 with the approval of the Japan Synchrotron Radiation Research Institute (JASRI) (proposal no. 2005B0561).