Mineral Resource Classification Using Geostatistical and Fractal Simulation in the Masjed Daghi Cu–Mo Porphyry Deposit, NW Iran

Abstract

Mineral resource classification is an important step in mineral exploration and mining engineering. In this study, copper and molybdenum resources were classified using a combination of the Turning Bands Simulation (TBSIM) and the Concentration–Volume (C–V) fractal model based on the Conditional Coefficient of Variation (CCV) for Cu realizations in the Masjed Daghi porphyry deposit, NW Iran. In this research, 100 scenarios for the local variability of copper were correspondingly simulated using the TBSIM and the CCVs were calculated for each realization. Furthermore, various populations for these CCVs were distinguished using C–V fractal modeling. The C–V log–log plots indicate a multifractal nature that shows a ring structure for the “Measured”, “Indicated”, and “Inferred” classes in this deposit. Then, the results obtained using this hybrid method were compared with the CCV–Tonnage graphs. Finally, the results obtained using the geostatistical and fractal simulation showed that the marginal parts of this deposit constitute inferred resources and need more information from exploration boreholes.

1. Introduction

Mineral resource classification is a significant step in mineral exploration projects for feasibility studies. At this stage, geostatistical methods can be used to support the assessment of uncertainty and risk. Resource classification can be carried out based on different criteria such as the distance between drillholes, kriging variance, multiple pass kriging plans, and uncertainty models that are provided by using geostatistical simulations [1,2,3,4,5,6,7,8,9].

Alternatively, one of the most important aims of geostatistical simulations is to assess the joint uncertainty between multiple realizations allowing a more complete representation of uncertainty for each block and between multiple locations of blocks, which deterministic approaches such as kriging are not able to produce. Geostatistical simulation methods reproduce input data parameters based on their statistical characteristics including their variogram and histogram [10,11,12,13,14,15,16,17,18].

Conventional geostatistical simulation approaches for continues variables (e.g., the ore grade) consist of the Sequential Gaussian Simulation (SGS) [19,20], the turning bands simulation (TBSIM) [21], the LU decomposition [22], the direct sequential simulation [23], the direct block simulation [24], and the probability field simulation [25]. The TBSIM was introduced by Matheron (1973) [26] as an unconditional simulation and developed by Journel (1974) [27] for the simplification of the Gaussian simulation problem in multi-dimensional spaces. This method produces a variable in one dimension and spreads them into 2D and 3D spaces for different realizations [28,29,30,31].

Fractal modeling is a proper methodology for the definition of various populations in natural sciences that is established by Mandelbrot (1983) [32]. The Concentration–Volume (C–V) fractal model proposed by Afzal et al. (2011) [33] has been used for the detection of different mineralized zones in different ore deposit types [34,35]. Furthermore, this method was used for the separation of ore grades’ populations from the different realizations [36].

In this research, the combination of the TBSIM and C–V model is utilized for copper and molybdenum resource classification in the Masjed Daghi Cu–Mo porphyry deposit, NW Iran. The coefficient of variations for Cu and Mo were computed based on the TBSIM and categorized by the C–V fractal modeling. This hybrid method between TBSIM and C–V fractal model based on the Conditional Coefficient of Variation (CCV) is a proposed method for this purpose. The average values of CCVs for Cu and Mo are used for resource classification. At last, the results obtained using this methodology are compared with the CCV–Tonnage graphs for Cu and Mo that were proposed by Dimitrakopoulos et al. (2009).

2. Geological Setting

The Masjed Daghi deposit is located 35 km east of Jolfa, East Azarbaidjan (NW Iran; Figure 1). This deposit is situated in the NW part of the Urumieh–Dokhtar magmatic belt (Arasbaran). This deposit contains 340 Mt with Cu and Mo averages equal to 2700 ppm and 60 ppm, respectively [37]. The host rocks include Eocene andesite, trachyandesite, dacite, tuffs, and agglomerates, intruded by a shallow diorite porphyry intrusion. Geological and mineralogical studies show that the quartz diorite subvolcanic unit along with Oligocene volcanic rocks in the area were severely altered by hydrothermal fluids. There are porphyry copper type and high sulfidation epithermal gold mineralization [37]. Six types of alteration zones were developed within magmatic rocks including potassic, potassic–phyllic, phyllic, phyllic–argillic, intermediate argillic, and propylitic [38].

Mineralization occurs in breccias and vuggy quartz and quartz–barite veins. Pyrite and chalcopyrite are the main ore minerals that associate with molybdenite, sphalerite, galena, and Fe-oxides. The gold contents vary between 0.1 and 30 ppm. The auriferous veins are surrounded by the argillic, silicic, alunite, and chlorite alteration zones. The epithermal system at the Masjed Daghi appears to be associated with a porphyry Cu–Mo system at deeper levels [38].

3. Material and Methods

3.1. Dataset

The dataset contains 18,679 core samples with 2 m length from 69 boreholes (Figure 2). These samples were analyzed using the ICP-MS method for Cu, Mo, and related elements at the ACME CHEMEX Company (Canada). The detection limits for Cu and Mo are 1 ppm and 0.5 ppm, respectively. The first step in exploratory data analysis includes the detection of all possible outliers and duplicated samples [41]. The outliers are recognized by the histograms of these data so there are no censored data and, finally, the outliers are corrected using the Durfel method [42].

3.2. TBSIM

The TBSIM was first introduced by Chentsov (1957) [43] in a special case of Brownian random functions. It has been developed to the stationary and stochastic simulation by Matheron (1973) [26]. This method aims at simplifying the simulation operation in multi-dimensional spaces according to simulation in one dimension and spreading them to the 2D/3D spaces [21]. The TBSIM is carried out based on a set of randomly distributed bands or lines. The general procedure follows [21,44]:

-

Raw data values should be declustered;

-

Declustered data are transformed to standard normal distribution;

-

Variograms are generated based on the transformed data;

-

Turning band simulation using the transformed data and variogram is performed;

-

Then, the simulation data must be back-transformed from normal distribution to the original scale;

-

Eventually, the realization maps of TBSIM are generated.

3.3. C–N Fractal Model

The concentration–number (C–N) fractal model was proposed by Hassanpour and Afzal (2013) [45] based on the number–size (N–S) fractal method [32]. The C–N fractal model is an essential fractal model used to define the geochemical anomaly or mineralized zone threshold values in an exploratory dataset. The C–N fractal model was defined as:

N(≥ρ) ∞ ρ^−β

(1)

where N(≥ρ) indicates the sample number with concentration values greater than the ρ. The ρ is the concentration of the element and β is a fractal dimension. The main advantage of this method is the classification of geochemical populations before estimation [46,47,48,49,50,51,52].

3.4. C–V Fractal Model

The C–V fractal model was first introduced by Afzal et al. (2011) [33] for the delineation of mineralized zones and barren host rocks in different types of ore deposits, especially porphyry ores. It must be added here that, in the C–V model, “C” can be replaced by either “concentration” (e.g., grade or tonnage) or “probability” (e.g., uncertainty). In this paper, the researchers used “C” to refer to both uncertainty and concentration. The C–V fractal model can be expressed as:

V(c ≤ ν) ∞ c ^−a1; V (c > ν) ∞ c ^−a2

(2)

where V(c ≤ ν) and V (c > ν) represent volumes (V) with concentration values (c) smaller and greater than contour values (ν), respectively; a1 and a2 are characteristic exponents [53,54,55]. In this paper, the Conditional Coefficient of Variation (CCV) values derived using a TBSIM are classified using the C–V fractal model.

3.5. Hybrid Regression Models

Different approaches have been proposed and developed for “Mineral Resource Classification” such as the number of drillholes and samples that are used in the estimation, probability density function of the grade for each sub-cell, and its estimated variance [3]. Dimitrakopoulos et al. (2009) [3] studied the estimated variances of a block model for a diamond deposit. In this method, the Conditional Coefficient of Variation (CCV) for each sub-cell is calculated using the following equation:

$$\mathrm{CCV}\left(\mathrm{u}\right)=\frac{\sqrt{{{\displaystyle \sum }}_{\mathrm{K}=1}^{\mathrm{K}+}{\left[{\overline{\mathrm{z}}}_{\mathrm{k}}-{\mathrm{z}}_{\mathrm{E}}^{*}\left(\mathrm{u}\right)\right]}^2\left[\mathrm{F}\left(\mathrm{u};{\mathrm{z}}_{\mathrm{k}}\right)-\mathrm{F}\left(\mathrm{u};{\mathrm{z}}_{\mathrm{k}-1}\right)\right]}}{{\mathrm{z}}_{\mathrm{E}}^{*}\left(\mathrm{u}\right)}$$

(3)

where k are threshold values for discretizing the range of variation of z values. z_k is the mean of the class z_k−1. In the case of a within-class linear interpolation model, it corresponds to the following:

z_k = (z_k−1 + z_k)/2

(4)

A first confidence interval based on the CCV involves fixed threshold CCVs derived using the fractal simulation that is compared to the CCV obtained using CCV versus the tonnage graph [3].

4. Results and Discussion

4.1. Pre-Processing

The first step of this methodology is using a dataset that is obtained using homotopic sampling patterns [56]. It is assumed that all conditions of the samples have equal length. Accordingly, the composites with lengths of 2, 5, and 10 m are created and the optimal length was selected as 10 m based on their histograms and comparison with the raw data. Based on this composition, the number of data was significantly reduced to 4057 samples.

In the next step, the data are declustered to ensure the global distribution is more representative using the cell declustering method [20,27]. The cell declustering technique is implemented in a dimension of 100 m × 100 m × 10 m based on the primary pattern of the drilling grid. The cell declustering technique [11] is applied in this study to correct the pseudo-skewness in the global distribution of Cu and Mo. The declustered data are then transformed to a standard normal distribution [20]. The histogram of raw, composited, and normal score data are depicted in Figure 3.

4.2. Application of TBSIM for CCV Calculation

The TBSIM approach in this study was applied to the CCV simulation of Cu and Mo in the Masjed Daghi deposit. To implement the TBSIM method, the variograms of Cu and Mo must be quantified over the transformed Gaussian values. The anisotropy studied by the experimental variogram in various directions, and the results indicated that the spatial distribution is isotropic; therefore, the omnidirectional variograms are calculated, and the two-structured model is fitted into the variables as depicted in

Table 1. Characteristics of fitted variograms for Cu and Mo.

Variable	Nugget	First Structure		Second Structure
		Sill	Range (Meters)	Sill	Range (Meters)
Cu	0	0.3	60	0.7	495
Mo	0.1	0.3	117	0.6	459

.

The TBSIM was carried out for 100 realizations with a MATLAB code that was provided by Emery and Lantuejoul (2006) [21] and the realizations were visualized using SGeMS software. The average (E-type) and variance of each realization for Cu and Mo indicate the proper reproduction of the mean and variance in the TBSIM as performed in Figure 4. Furthermore, the variograms of different realizations were compared with the variograms of the raw data for the validation of the TBSIM operation as depicted in Figure 5. The CCVs were calculated and visualized for each realization of Cu and Mo (Figure 6). Finally, the average of these CCVs was calculated for Cu and Mo as E-types (average realizations), which are represented in Figure 7.

4.3. Mineral Resource Classification using CCV–V Fractal Modeling

First of all, the C–N fractal modeling was carried out on Cu and Mo raw data for the detection of the main thresholds (Figure 8). These main thresholds were used for the separation of the Cu and Mo mineralized zone from barren host rocks. These thresholds were calculated 0.25% and 60 ppm for Cu and Mo, respectively. The sub-cells with Cu and Mo grades lower than these thresholds were removed from the resource classification as barren host rocks. On the other hand, the mineral resource classification was carried out for sub-cells with Cu ≥ 0.25% and Mo ≥ 60 ppm as the main mineralized zone.

Based on the calculated CCVs in different sub-cells, the CCV–V log–log plots were generated for the Cu and Mo CCVs E-types, as depicted in Figure 9. The different populations were separated according to the changing dips of the segments [35,49]. There are three and four populations for CCVs of Cu and Mo, respectively. Based on these log–log plots, the Measured, Indicated, and Inferred copper resources have CCVs < 0.18, 0.18–0.35, and >0.35, respectively (

Table 2. Mineral resource classification by CCV–V fractal modeling for Cu and Mo.

Resources Type	CCV Range for Cu	CCV Range for Mo
Measured	<0.18	<0.84
Indicated	0.18–0.35	0.84–0.96
Inferred	>0.35	>0.96

). Moreover, there are Measured, Indicated, and Inferred Mo resources that contain CCV values equal to <0.84, 0.84–0.96. and >0.96, respectively. On the other hand, the first and second populations in these log–log plots were defined as Measured and Indicated resources (Figure 9 and Table 2).

Based on study by Dimitrakopolous et al. (2009) [3], the CCV–Tonnage graphs for Cu and Mo were generated, as depicted in Figure 10 and

Table 3. Mineral resource classification using CCV–Tonnage graphs for Cu and Mo.

Resources Type	CCV Range for Cu	CCV Range for Mo
Measured	<0.2	<0.9
Indicated	0.2–0.4	0.9–1
Inferred	>0.4	>1

. These graphs were generated based on different realization of Cu and Mo (Figure 10). The resource classification was carried out based on the turning points in the CCV–Tonnage plots for Cu and Mo. The Measured and Indicated resources for Cu and Mo are larger than the results derived using fractal modeling. Their CCV values for the Measured and Indicated resources are <0.2 and 0.2–0.4 for Cu, respectively (Table 3). However, the CCVs variations for Mo resources are <0.9 and 0.9–1 for Measured and Indicated classes, respectively (Table 3). The location of different classes were shown in Figure 11 and Figure 12. The classes resulting from using fractal modeling have better correlations with the density of the boreholes.

5. Conclusions

In this study, the TBSIM and CCV–V fractal model were used for the delineation of different mineral resource classes in the Masjed Daghi Cu–Mo porphyry deposit (NW Iran). The results obtained using this methodology show that there are smaller Measured and Indicated classes derived using the CCV–V fractal model in comparison with the results of the CCV–Tonnage graphs. The measured resources for Cu are <0.18 and <0.2 based on the hybrid and CCV–Tonnage graph, respectively. Furthermore, this resource class is <0.84 and 0.9 as derived using the hybrid and CCV–Tonnage graph, respectively. Consequently, many sub-cells in the marginal part of this deposit are classified in the Indicated level using the CCV–Tonnage graph but these sub-cells are categorized to the Inferred class based on the hybrid of TBSIM and CCV–V fractal model. There are low values of borehole data, which reveals the higher accuracy of the hybrid method. On the other hand, the hybrid methodology of the geostatistical simulation and CCV–V fractal modeling can be used for mineral resource classification in different types of ore deposits. The CCV–V fractal modeling can be separated into different resource classes with higher accuracy in comparison with other methods.