A rapid and precise quantitative electron probe chemical mapping technique and its application to an ultrahigh-pressure eclogite from the Moldanubian Zone of the Bohemian Massif (Nové Dvory, Czech Republic)

An electron probe microanalyzer (EPMA) is a common tool to analyze chemical compositions of minerals. Its use in modern petrology allows for quantification of chemical heterogeneities or textures within rocks, thereby enabling the rock history to be more accurately determined. By describing chemical textures of rocks at an appropriate scale, their pressure–temperature evolution and history of mass transfer may be deduced.

Two common applications of an EPMA are quantitative spot analysis and X-ray mapping. The former collects an X-ray signal from an unknown sample that is compared to a standard signal intensity to quantify elemental concentrations. The latter is employed to create X-ray intensity images (X-ray maps) to describe the qualitative distribution of elements. X-ray maps also show relationships between micro-textures and variations in the chemical composition of minerals (zoning), which reflect the conditions of crystallization (e.g., Marmo et al. 2002; De Andrade et al. 2006; Vidal et al. 2016; Ganne et al. 2012; Lanari et al. 2014; Trincal et al. 2015; Lanari and Engi 2017). Theoretically, X-ray intensities collected from a sample surface are proportional to the concentration of the elements in an excitation volume for both point and mapping analyses. As such, quantitative X-ray mapping is a possible application of EPMA analysis. However, creating a quantitative X-ray map by gathering quantitative point analyses requires an unrealistically long period of time because a typical quantitative point analysis requires ca. 1 min for 10 elements (e.g., Si, Ti, Al, Cr, Fe, Mn, Mg, Ca, Na, and K). Under these conditions, approximately one week would be needed to create a 100 × 100 pixel map.

Previous studies have proposed other means to construct quantitative X-ray maps within a reasonable time frame (Tracy et al. 1976; Clarke et al. 2001; Tinkham and Ghent 2005; De Andrade et al. 2006; Lanari et al. 2014, 2018). Among them, XMapTools (Lanari et al. 2014, 2018) provides a reliable analytical procedure by combining mapping results obtained at a higher current intensity with referential spot analysis measured within a mapped area. XMapTools addresses matrix effects by implementing (1) classifying mineral phases based on mapping data and (2) converting X-ray intensities into calibrated mass concentrations for each mineral phase by comparing the X-ray map to referential spot analysis within a mapped area. Both steps assume each pixel of the mapped area is attributed to a single mineral phase and require a sufficiently small mapping-beam diameter. Although XMapTools is powerful in the analysis of micro-textures or coarse-grained centimeter-scale textures (Lanari and Engi 2017), a problem arises when the X-ray map employs a large step size to reduce mapping time. When a mapped area contains fine-grained textures such as accessory minerals and symplectite, a small beam diameter with a large step size may miss grains that were analyzed as an internal standard, particularly for inequigranular textures. In addition, such low-resolution maps are inadequate for use in discussing variations in local bulk compositions because of possible underestimation of the abundance of fine-grained phases.

This study provides a method to quantify X-ray maps with a scalable probe diameter, which is particularly powerful in analyzing centimeter-scale areas containing fine grains within a reasonable time frame (within a few days in the present case study but potentially within a day). Enlarging the probe diameter increases the number of pixels that traverse the grain boundaries (i.e., multi-phase pixels) that cannot be treated properly under the assumption that each pixel is to be attributed to a single phase. To overcome this limitation, we introduced a distribution-based cluster analysis (Titterington et al. 1985), an unsupervised machine learning technique, to evaluate the mixed content in such pixels. Evaluating the content of the multi-phase pixels enables reliable estimation of the relationship between the X-ray map and referential spot analysis. The following sections describe the principle of the new method to compile reliable quantitative maps. As a case study, we applied the proposed method to an ultrahigh-pressure eclogite collected from the Moldanubian Zone of the Bohemian Massif (Nové Dvory, Czech Republic). The results show that centimeter-scale quantitative chemical mapping is a powerful tool for evaluating millimeter-scale compositional variations in rocks.

The fundamental idea behind our quantitative electron probe chemical mapping technique follows internal standardization methods (De Andrade et al. 2006; Lanari et al. 2014, 2018). Thus, our technique requires referential spot analysis of all mineral phases within an area to be mapped, X-ray mapping of the area (both using EPMA), and additional data processing by computer. During the data processing, a cluster analysis was performed on map data to address variation in the matrix effects, after which a regression analysis was used to convert mapping X-ray intensities to mass concentrations. The data processing was completed using QntMap, an open source R package developed in this study and distributed on a social coding platform, GitHub (https://github.com/atusy/qntmap accessed on February 10, 2018). QntMap was designed to address text files in the format of the output by JEOL EPMA. However, the proposed algorithm can be applied to any EPMA results by suitably preparing the required files (Supplemental¹ Table S1).

This study employed an EPMA with five wavelength-dispersive X-ray detectors (JEOL JXA-8105) at the Department of Geology and Mineralogy, Kyoto University, Kyoto, Japan. Quantitative point analysis was performed with 15 kV of accelerating voltage, 10 nA of probe current, 3 μm of probe diameter, and 10 and 5 s of dwell time for counting the peak and background, all of which are ordinary employed analytical conditions. X-ray mapping was performed with 15 kV of accelerating voltage, 100 nA of probe current, 15 μm of probe diameter, 20 μm of step size, and 120 ms of dwell time. The probe current and the dwell time follows that of De Andrade et al. (2006) and enables high precision without significant damage of matrices. The probe diameter was slightly smaller than the step size to account for the X-ray excitation volume. Dead time corrections for mapping data were performed using software in EPMA.

X-ray maps inevitably contain a certain number of multi-phase pixels, whose chemical compositions are a combination of chemical compositions of the phases weighted by their abundance. In our method, the abundance of phases in each pixel is estimated by a soft cluster analysis (Fig. 1), which assigns membership grades to each data point to express ratios of clusters. The membership grades are interpreted in terms of mass ratios of the phases in the multi-phase pixels because the cluster analysis is applied to mapping data reflecting chemical compositions (X-rays and backscattered electrons). This approach is in contrast to hard cluster analyses, such as the k-means clustering algorithm (KCA; MacQueen 1967) employed by XMapTools, which assigns exactly one cluster to each data point. Although a KCA is a simple and commonly used algorithm and is known to be a powerful tool for geochemical studies (e.g., Iwamori et al. 2017), it results in overestimation or underestimation of chemical components in multi-phase pixels. For example, consider a multi-phase pixel consisting of 75% quartz (SiO₂) and 25% rutile (TiO₂) as shown in Figure 1. Hard cluster analysis assigns such pixels to a single cluster of quartz and provides 100% SiO₂ as the corresponding composition. In contrast, a soft cluster analysis evaluates the membership of the constitutive clusters and allows classification of the pixel as a combination of mixed phases.

Among various soft-clustering algorithms, QntMap implements a distribution-based clustering algorithm (DBCA) (Titterington et al. 1985). In a DBCA, a data set is considered to comprise a given type of multivariate mixture probability distribution (e.g., Gaussian or Poisson) and clusters are considered to be composed of groups of data belonging to each distribution in the mixture. Estimations of distribution parameters and data point assignment to clusters are performed iteratively until convergence is achieved. QntMap employs the implementation of DBCA described in Witten (2011) that assumes input data comprise a mixed Poisson distribution, which is appropriate for EPMA maps that count photons of X-rays and backscattered electrons.

As an initial guess, data from a referential spot analysis are set to initial cluster centers, and mapping data points are assumed to belong to the nearest cluster centers. This assumption enables rapid convergence over a small number of iterations. The DBCA algorithm we employed (Witten 2011) includes a power transformation to account for overdispersed data. Thus, solid-solution phases with small compositional variation would be correctly partitioned into similar clusters. However, phases with significant compositional variation or zoning (e.g., >5 wt%) or multi-modal patterns on a histogram may be treated as multiple phases, such as garnet1 and garnet2.

In addition to evaluating membership degrees, a DBCA overcomes other limitations in a KCA. A KCA is a special case of a DBCA that assumes that a data set is composed of a mixture of multivariate Gaussian distributions that have an isotropic shape in multidimensional space, do not overlap, and whose integrals over the entire space yield similar values. Thus, a KCA does not work well when a population of clusters varies significantly (Guha et al. 1998), as shown in Figure 2a, which is a common case for modal compositions of rocks. The use of a KCA is also inappropriate when the shapes of groups are anisotropic (Nagy 1968), as shown in Figure 2b, which may arise because of solid solutions in minerals or differences in the sensitivities of the X-ray detectors. These factors commonly result in a KCA partitioning major phases into too many clusters and minor phases into clusters representing different mineral phases. A DBCA is generally free of the assumptions made in a KCA (Fig. 2) and is better able to distinguish mineral phases (Fig. 3).

Quantitative EPMA analyses are based mainly on a linear relationship between the X-ray intensities of unknown samples and known standards, employing well-established empirical models such as ZAF (e.g., Heinrich 1991). The relationship between concentrations and X-ray intensities for the i-th element is usually written as follows:

where, C represents the mass concentration, I represents the X-ray intensity, and G is the correction factor for the atomic number, absorbance, and fluorescence derived from coexisting elements. The superscripts “unk” and “std” indicate unknown samples and standard materials, respectively. The subscripts “i,” “net,” and “qnt” indicate the element, the net X-ray intensity (calculated as peak minus background), and quantitative analysis, respectively. Previous studies have showed that this simple linear relationship can be also adapted to a combination of X-ray maps and referential spot analyses performed within a mapped area (De Andrade et al. 2006; Lanari et al. 2014, 2018). These studies converted mapped peak X-ray intensities to concentrations by using the following calibration curve:

where, A is a correction factor based on the mineral phases and elements analyzed. The subscripts “pk,” “bg,” and “map” indicate the peak X-ray intensity, the background X-ray intensity, and mapping analysis, respectively. The term A represents the ratio of the results of the referential spot analyses (⁠$Ciunk$⁠) in Equation 1 and the net X-ray intensity of the mapping analysis; thus, the effects of the atomic number correction and other factors (i.e., G in Eq. 1) are included.

Previous studies have proposed background corrections. De Andrade et al. (2006) directly analyzed background intensity with mapping conditions; however, this approach requires additional mapping analysis. Donovan et al. (2016) employed a mean atomic number (MAN) correction that required X-ray mapping of not only the elements of interest but of all the major elements. XMapTools (Lanari et al. 2018) used the variability observed in mass concentrations of internal standards, though it cannot be applied to phases/elements of low chemical variability, i.e., minor elements such as chromium in mafic rocks.

Considering the effects of multi-phase pixels and background corrections, the fundamental idea of Equation 2 is as follows:

where, X is the abundance of a phase in a pixel whose sum is 1, α is the ratio of mass concentration and net X-ray intensity from spot analysis, β is a correction factor relating the difference in analytical conditions between spot and map analyses, and γ is a correction factor for background X-ray intensities. Refer to Supplemental Material A¹ for the derivation of Equation 3. As the membership degrees are interpreted as mass ratios, not volume ratios, density correction was considered to be discarded.

In Equation 3, α is a correction factor dependent on the mineral phase and element and represents the ratio of mass concentration to quantitative net X-ray intensity (i.e., the slope of the regression lines shown in Fig. 4). In other words, α corresponds to $Cistd$G / $Ii,net,qntstd$ in Equation 1. Although α depends on matrices, it can often be approximated as a constant for each phase as shown in Figure 4, in which the relationship between the net X-ray intensities of the referential spot analysis and the corresponding result of MgO concentrations are depicted. For each phase, the relationship can be explained by a single straight line passing the origin. When phases show considerable variation in MgO concentration, MgO-poorer data points tend to occur on the left side of the regression lines and MgO-richer points on the right side (Figs. 4c–4d). Yet, the regression lines overlap with the 95% confidence intervals of each data point, indicating the relationships can be explained by single calibration curves for each phase.

In Equation 3, β is a correction factor that depends on the element, representing the ratio of quantitative to qualitative peak X-ray intensities (i.e., the slope of the regression line shown in Fig. 5). Its product with α (= αβ) corresponds to A in Equation 2. β is independent from the matrix effects because the ratio of quantitative to qualitative peak X-ray intensity is determined from the difference in probe current, dwell time, and probe diameter. The difference in probe diameter is less significant than that of the others, however, as a larger probe diameter decreases X-ray intensities because of the geometrical relationships between the X-ray excited region and X-ray detector. Figure 5 shows the relationship between MgO peak X-ray intensities from referential spot analysis and X-ray mapping from identical coordinates within a natural sample. Most data points occur on a single line, regardless of mineral species. The outliers are fine-grained phases within symplectites, such as amphibole, diopside, and plagioclase. Mapped pixels among the outliers are analyzed multiple phases in addition to the that targeted by referential spot analysis. As outliers belong to multi-phase pixels, they yield low membership grades of clusters that correspond to phases analyzed in the referential spot analyses. Thus, β is estimated for each element by a weighted least-squares regression that weighs data points by the membership grade of the phase analyzed in the referential spot analyses, or by 0 for those that are regarded as very fine grains during the exploration of spots analyzed as a reference. Proper weighing enables confident estimation of β.

The term γ is a correction factor for background X-ray intensities, estimated as

and it depends on both the element and phase. N is the number of quantified points for a certain phase and j indicates the j-th quantified point. This formulation assumes that the matrix effects on backgrounds are approximately constant. γ values divided by β values (= γ/β) correspond to $Ii,bg,mapphase$ in Equation 2. Figure 6 shows the relationships between concentrations and background intensities for MgO in various phases. Although least-squares regression curves (blue dashed lines) show positive or negative slopes and indicate the presence of matrix effects on backgrounds, the average value indicated by the red line overlaps with the 95% confidence intervals of most points. Thus, we employed average values for the background of each phase as γ. Our approach provides similar background values to those of Lanari et al. (2018) when considerable variations are recognized in the mass concentration of the elements of interest.

As $βiIi,pk,mapunk,phase$ and $γiphase$ are equivalent to the peak and background X-ray intensities under quantitative analytical conditions, respectively, $βiIi,pk,mapunk,phase—γiphase$ provides the net X-ray intensity under X-ray mapping conditions. This formulation is considered to be valid because both peak and background X-ray intensities are proportional to the probe current. In addition, this formulation is robust against multi-phase pixels. As shown in Figure 7, the calibration curves of omphacite (a major phase of the studied sample) and diopside (mostly observed as fine-grained symplectite) are based on both Equations 2 and 3. They are expected to be identical because both omphacite and diopside belong to the clinopyroxene group and the relationship between mass concentrations and mapping peak intensities is fairly similar. Our proposed method provides approximately identical curves (Fig. 7b), while the simple least-square calibration using Equation 2 provides different calibration curves for diopside and omphacite because of the abundant multi-phase pixels containing diopside and other Mg-poor phases (plagioclase).

Our formulation estimated β, which is the standard-to-mapping intensity ratio, can be affected by multi-phase pixels including the abundance of phases in a pixel. Therefore, β for the multi-phase pixels was correctly estimated in the proposed method. β for certain elements can be estimated if there is any coarse-grained phase that contains the element. For example, β_MgO can be estimated even if fine-grained phases such as amphibole, diopside, and plagioclase are neglected (the black line shown in Fig. 5). This means that the proposed method can evaluate the compositions of fine-grained phases such as accessory and symplectite minerals even if the coarse grains that are measured as single-phase pixels are not found in the sample.

Here, we provide an example of an application of QntMap by analyzing a millimeter-scale bulk chemical variation in an ultrahigh-pressure eclogite collected from Nové Dvory in the Moldanubian Zone of the Bohemian Massif, Czech Republic. Nové Dvory is within a former continental collision zone and comprises ultrahigh-pressure metamorphic rocks (Medaris et al. 2005). Eclogite lenses and garnet pyroxenite layers occur within a block of garnet peridotite that is surrounded by felsic gneisses (Fig. 8). The eclogite outcrop consists mainly of a garnet-rich matrix and subsequent pyroxene-rich layers with a thickness of a few millimeters to a couple of centimeters (Medaris et al. 1995).

The study sample is a bimineralic eclogite that contains a 3 mm thick pyroxene-rich layer. Garnet, omphacite, and rutile represent an eclogite-facies mineral assemblage in the sample (Fig. 3). The garnet-rich matrix contains ca. 20–70 vol% garnet and 10–50 vol% pyroxene (i.e., omphacite + diopside). The pyroxene-rich layer contains up to 10 vol% garnet and 60–70 vol% pyroxene. Garnet is partially replaced by plagioclase, amphibole, diopside, spinel, and hematite, and omphacite is replaced by diopside and plagioclase. Rutile is commonly transformed to ilmenite. The degree of retrogressive replacement is more extreme in the pyroxene-rich layer, in which eclogite facies minerals comprise only 30–50 vol%, compared to 50–80 vol% in the garnet-rich matrix (Fig. 9). Pyroxene and plagioclase are intercalated in the pyroxene-rich layer and are finer-grained compared to the garnet-rich matrix.

We applied QntMap to an analyzed area of 33 × 10 mm in a polished thin section of the studied sample (Fig. 3). The EPMA mapping analysis took 55 h. The validity of the calculated quantitative chemical maps was strongly dependent on the reliability of the calibration curves, which were controlled by the number of referential spot analyses. XMapTools empirically requires a minimum of 20 referential spot analyses for each phase present (Lanari et al. 2014). In this study, the required number of referential spot analyses was statistically determined by a Monte Carlo (MC) cross-validation (CV; Supplemental¹ Fig. S4). First, the sample data were divided into a training data set and a test data set. Then, a calibration curve was calculated using the training data set. We evaluated the validity of the calculation by assessing the degree to which the calibration curve satisfied the variation in the test data set. The validity of this approach was assessed using the root mean square error (RMSE) calculated between test data and estimated values. A small RMSE indicates that the calibration curve reasonably explains the relationship among variables. In MCCV, CV iterations are performed using a randomly chosen number of training data. We performed MCCV with up to 10 000 iterations for each size of training data set to estimate α and β in Equation 3. Figure 10 shows the result of MCCV for MgO, using a total of 159 and 354 point analyses for the α and β values, respectively. For both α and β, the variance and the mean RMSE decrease and approach a constant value with increasing size of the training data set at approximately 15 analyses. Therefore, we concluded that a minimum of 20 point analyses of each phase is adequate for analysis using QntMap.

Based on the quantitative chemical maps obtained by Qnt-Map, a compositional profile across the layering structure is estimated by averaging the compositions of pixels in the same column, weighted by the density of phases in each pixel (Fig. 9). Approaching the pyroxene-rich layer, the garnet-rich matrix shows increases in Cr₂O₃ (0.0–0.3 wt%), MgO (9–14 wt%), and X_Mg [Mg/(Fe+Mg) = 0.5 to 0.8] and decreases in FeO (15–5 wt%). Note that these chemical gradients occur over a shorter distance for Cr₂O₃ (∼2 mm) than for MgO and FeO (∼12 mm). Concordantly, the Cr₂O₃ content of garnet and omphacite increases from the garnet rich layer (<0.5 wt% in both garnet and omphacite) toward the pyroxene-rich layers (up to 4 wt% in garnet and 2 wt% in omphacite) according to spot analysis.

Figure 11 shows the X_Mg and Cr₂O₃ contents of (local) bulk compositions from the studied sample and eclogite and pyroxenite from the Gföhl Unit, a geological unit of Nové Dvory (Beard et al. 1992; Medaris et al. 1995; Obata et al. 2006), and gabbroic rocks from the South Indian Ridge (Niu et al. 2002). Representative local bulk compositions of the garnet-rich matrix and pyroxene-rich layer in the Nové Dvory sample were also estimated from the shaded areas shown in Figure 9 (Table 1) and were plotted as shown in Figure 11. These were chosen to avoid chemical gradients. QntMap shows that the garnet-rich matrix corresponds well to the Gföhl eclogite and the gabbroic rocks, whereas the pyroxene-rich layer corresponds to the Gföhl pyroxenites. Medaris et al. (1995) distinguished pyroxenite from eclogite on the basis of elemental compositions, with X_Mg > 0.80, Na₂O < 0.75 wt%, Cr₂O₃ > 0.15 wt%, and Ni > 400 ppm. Based on these criteria, the garnet-rich matrix of the study sample corresponds to the eclogite (X_Mg = 0.61, Na₂O = 1.4 wt%, and Cr₂O₃ = 0.03) and the pyroxene-rich layer is intermediate in character between the eclogite and pyroxenite (X_Mg = 0.79, Na₂O = 1.60 wt%, and Cr₂O₃ = 0.26 wt%).

Previous studies have suggested that the Nové Dvory eclogite accumulated from a high-pressure melt (Medaris et al. 1995) or was transformed from a lower-pressure gabbroic cumulate (Nakamura et al. 2004; Obata et al. 2006). However, the chemical compositions of garnet and omphacite differ significantly between the garnet-rich matrix and pyroxene-rich layer. Perhaps the layering structure was produced by metasomatism triggered by different source materials, and thus we suggest a combination of the aforementioned two hypotheses; i.e., that the Nové Dvory eclogite was transformed from a subducted gabbroic rock that was partially intruded and metasomatized by pyroxenitic melt under high-pressure conditions. This scenario explains the rare occurrence of pyroxene-rich layers in the Nové Dvory eclogite and the millimeter-scale variations in X_Mg and Cr₂O₃. The chemical gradient occurred as a result of diffusion and corresponding chemical reactions during intrusion of pyroxenitic melt. The chemical gradient of MgO occurred along a longer distance than that of Cr₂O₃ because of the higher diffusivity of Mg compared to Cr and possible later modification by fluid activities during amphibolitization. Previous studies identified millimeter- to centimeter-thick pyroxenite layers from the host peridotite body and concluded that they were high-pressure cumulates from a basaltic melt based on trace elements and isotopes (Medaris et al. 1995; Svojtka et al. 2016). These inferences are also consistent with our data.

QntMap potentially provides high-quality mineral distribution maps and quantitative chemical maps at a centimeter cale based on EPMA data. Scalable lateral resolution in X-ray maps enables an increase in the probe current, as less lateral resolution produces less damage to the sample at a fixed probe current and dwell time. This means that the rate of qualitative mapping can be further increased from the current application (18 h/cm²) by reducing the dwell times at higher probe currents. QntMap is also powerful in millimeter-scale mapping of fine-grained textures such as symplectites, because multi-phase pixels occur when lateral resolutions are not sufficiently high. The maximal lateral resolution is defined by the quantitative point analyses, typically a few micrometers. For bulk chemistry, X-ray fluorescence (XRF) is a conventional technique; however, it is a destructive method and is not easily applied to millimeter-scale compositional domains within rocks. QntMap preserves textural information and thus is powerful not only when domains are very small, but also when domains are transitioning or showing complex shapes. QntMap is a highly scalable technique to create quantitative electron probe chemical maps.

The application of QntMap was successful in quantitatively evaluating continuous changes in local bulk compositions and textures at a millimeter scale that tends to be overlooked in bulk-rock chemistries of hand-specimen samples. The present results suggest that the layering structure in the Nové Dvory eclogite was produced by intrusion of a pyroxenitic melt into a subducted gabbroic rock at high-pressure conditions. The preservation of millimeter-scale inhomogeneities in the Nové Dvory eclogite implies the duration of high-temperature conditions (>1000 °C in peak) was sufficiently short to avoid complete homogenization by volume diffusion. This confirms rapid exhumation of the eclogite (Nakamura et al. 2004) and its host peridotite (Medaris et al. 1990) as inferred from garnet zoning patterns, and it further implies the presence of a low-viscosity pyroxenitic melt that may have aided in the rapid exhumation of the metamorphic rocks. The formation process of millimeter-scale inhomogeneities may record a key event in petrogenesis and can be retrieved by combining millimeter-scale local bulk chemistry with textural information based on quantitative EPMA mapping.

This paper benefited from the careful and productive reviews of P. Lanari and J.W. Singer and editing by T. Mueller. This work was partly supported by the Inter-Graduate School Program for Sustainable Development and Survivable Societies of Kyoto University to A.Y., and JSPS KAKENHI Grants Nos. JP25257208 to T.H., JP16K17835 to K.Y., and JP25120005 and JP15K20864 to T.K. Financial support was also provided by JST PRESTO (Grant No. JPMJPR1676) to T.K. The Scientific Program RVO67985831 of the Institute of Geology, The Czech Academy of Sciences, to M.S. is also acknowledged.