Geostatistical Analysis for delineating sterile inclusions in Sidi Chennane ' phosphatic series , Morrocco

RESUMEN

The anomalies detected in phosphatic series of Sidi Chennane, one of phosphatic basins in Morocco, hinder the proper exploitation of phosphate levels and the assessing phosphate reserves seems incorrect.The purpose of this study was the evaluation of geostatistical methods for delimitation of these disturbances.To cover all the zones being able to be disturbed, we carried out, during the geophysical prospection in a parcel of 50 ha, 5151 resistivity measurements as horizontal profiling using the well-known Schlumberger array, in order to map the spatial distribution of the sterile hardpan inclusions.Geostatistical tools were used to quantify the spatial correlation between apparent resistivity data.Semivariograms were obtained using a classical Matheron semivariogram estimator and fit to the experimental semivariograms obtained.We have selected those with the best fit in terms of sum of squared residuals (SSR).Geostatistical analysis was performed using the software VESPER 1.63.Spatial distribution maps were made by ordinary kriging, the qualitative interpretation of these maps reflects that the exponential model is found to be the best model representing the spatial variability of our geoelectric data.The qualitative interpretation of the kriged resistivity maps allows defining resistivity contrast, consequently we have delimited the crossing dominate area from a "normal" into a "disturbed" area.Models of the geology were successfully obtained from geostatical method, which help mapping the phosphate deposit inclusions and the estimations of phosphate reserves were improved and better constrained.

Introduction
Morocco is the world's third largest phosphate producer, after the USA and China.Total mine production recorded by the Ministry of Energy and Mines in 2013 was 30 Mt so more than 75% of world reserves.Four major phosphate basins are now known and are being exploited, three of which are located in central-northern Morocco.The four main deposits of phosphate are: the Oued Eddahab basin situated in Sahara, the central ganntour basin near Youssoufia, the Meskala basin at east of Essaouira and the Oulad Abdoun basin situated near Khouribga.The existence of morocco sedimentary phosphate rock has been known since 1908 in the Meskala basin, but it had not generated significant interest until the discovery, in 1917, of the Oulad Abdoun basin.
The geological investigations carried out in Sidi Chennane phosphatic deposit in the Oulad Abdoun basin revealed a phenomenon of sterile hardpan [Bakkali 2005, Baba & al. 2012;2014a and2014b] these sterile bodies are formed by accumulations of silicified limestones or by limestone blocks within an argillaceous matrix which interfere with phosphate extraction and their resistivity, above 200 Ohm.m, is higher than the phosphate-rich mineral resistivity wich is near 120 Ohm.m [Bakkali and Amrani 2008;Baba & al. 2012].The application of the electric prospection methods constitutes a suitable means to map these sterile bodies [Baba & al. 2012[Baba & al. , 2014a[Baba & al. and 2014b]].
One of tools for mapping of electrical resistivity data is geostatistical methods.Resistivity data were analyzed using geostatistics (variogram) in order to describe and quantify the spatial continuity.
This study has been carried out with the aim of evaluation geostatistical methods for estimation of subsurface layers delimitation and specified their spatial variability in order to establish a model of sterile bodies' distribution and would permit the definition of these structures before the mining front reaches them.

Data acquisition
Electrical and electromagnetic (EM) methods have been important in the field of Applied Geophysics for about a century, particularly for shallow and near-surface investigations.Geoelectrical resistivity surveying has become an important and useful tool in hydrogeological studies, mineral prospecting and mining, as well as in environmental and engineering applications [Griffiths & al. 1990;Griffiths and Barker, 1993;Dahlin and Loke, 1998;Olayinka, 1999;Olayinka and Yaramanci, 1999;Amidu and Olayinka, 2006;Baba & al. 2012;Ouadif & al. 2012].
Electrical sounding is a method to investigate the change in earth resistivity with depth at a particular location.Traditionally, arrangements using four electrodes (two current-transmitting electrodes and two voltagesensing electrodes) are used for either vertical soundings or horizontal profiling.For vertical soundings, the electrodes are arranged symmetrically according to a center, with increasing distances between electrodes used to explore deeper depths [Baba & al. 2012;Ouadif &al. 2012 and2014a].
When using geoelectrical surveying techniques, an estimate of the earth resistivity is calculated by using the well-known relation between resistivity, an electric field, current density (called Ohm's Law), the geometry and spacing of the current and potential electrodes.When the earth is not homogeneous neither isotropic, this estimated value is called the apparent resistivity, p_app, which is an average of the true resistivity in the measured section of the earth.To determine the intrinsic resistivity for the different terms of the phosphatic series, we carried out, during the geophysical prospection in a parcel of 50 ha,5151 resistivity measurements as horizontal profiling (SYSCAL-R2 resistivity instrument IRIS Instruments) using the wellknown Schlumberger array, in order to map the spatial distribution of the sterile hardpan inclusions.The 5151 stations of the resistivity survey is a compilation of 51 profiles spaced at 20 m.For every profile, there were 101 stations separated by 5 m using the well-known Schlumberger array.

Geostatistics
Geostatistics, started with the work in geology and mining of Krige (1951) then developed by Matheron (1963) with his theory of regionalized variables, are methods of interpolation which predict unknown values from data observed at known locations, and it minimizes the error of predicted values which are estimated by spatial distribution of the predicted values.As variables are supposed having an effect on each other, the first concept is related to the variogram of values we want to estimate.It relates the similarity or difference, expressed as the semi-variance, between values at different places to their separation distance (lag) and direction.Distance that variables have an effect on each other is called range of a variogram and is an important concept because after that distance, variables have no effect on each other.
The plot of the semivariance versus the lag h is called semivariogram.In some cases, it increases with h and reaches an asymptotic value, a sill.Range is the value of h at which the semivariogram asymptotes.Values of the properties at sites separated by a distance shorter than the range will be spatially correlated.The range depends on the scale of the study and on the property evaluated.[Krige, 1951;Webster and Oliver, 2001].Let be r(x) and r(x+h) the values of electrical resistivity sampled at.
positions x and x+h, with h a distance, known as lag, which separates both positions, semi-variances can be estimated as: (1) For discrete sampling sites, such as earthen material samples, the function is estimated as: (2) Where r(x i ) : the value of electrical resistivityr at location of x i N(h) : the number of point pairs for distance class h.

Power: (6)
With h: lag distance c: sill a: range Nugget effect is essential to optimize the process under study: the higher its difference for the semivariogram sill, the higher the phenomenon continuity, the lower the estimation variance, or the higher the estimation confidence.
The variogram plot is fitted with a theoretical model which provides information about the spatial structure as well as the input parameters for kriging interpolation.Kriging algorithms (Krige, 1951) are described as a best linear unbiased estimator (BLUE), which is a process of a theoretical weighted moving average: (7) where r(x 0 ) is the value to be estimated at the location of x 0 , r(x i ) is the known value at the sampling site xi.And l i denotes the weight of the i th observation, they are calculated searching an unbiased estimator of z with minimum variance [Webster, 1985].
The estimation errors (or kriging variances) need to be minimised.To ensure that the estimate is unbiased, the weights need to sum to one:

Statistical Analysis
The classical statistical procedures assume that variations are randomly distributed within each class of observations and independent of sample location [Mohmmadi, 2006].The mean is therefore used for the estimation of properties for unsampled locations within each class.Statistics of dispersion (e.g.variance, standard deviation, coefficient of variability, confidence limits) are used to show the precision of the mean as an estimator.The resistivity distributions are not normal, which is the particularity of the anomalous bearing population.Variation between maximum and minimum for these data shows a wide range.Assuming the samples were collected randomly from different locations in the field, we can conclude from these statistics that resistivity at the upper sublayers is more variable on the field than that at the lower layers.

Geostatistical methods
Vesper 1.6 software (ACPA, University of Sydney) developed by the Australian Center for Precision Farming [Minasny & al. 2002] was used to calculate experimental electrical resistivity semivariograms of experimental fields in order to obtain the parameters: range, sill and nugget variance for evaluating the spatial correlation of these variables.Semivariograms were obtained using a classical Matheron semivariogram estimator.
The experimental semivariograms were fitted with various theoretical models like spherical, exponential, gaussian, linear and power by an automatic (least-squares) estimation of semivariogram parameters.Indicator variograms were built adopting a lag of 50 m.The theoretical model that gave minimum standard error is chosen for further analysis.
Table 2 shows semivariogram models and parameters.The attributes were fitted to the exponential model.The ratio of nugget to total semivariance, expressed as a percentage, was used to classify spatial dependence [Cambardella & al. 1994], a ratio less than 25% indicated strong spatial dependence, between 25 and 75% indicated moderate spatial dependence and more than 75% indicated weak spatial dependence.Spatial interpolation accuracy and precision were evaluated through Spatial interpolation accuracy and precision were evaluated through crossvalidation approach.Most important criteria include: mean bias error (MBE), mean absolute error (MAE), mean squared error (MSE) and root mean squared error (RMSE).MAE and RMSE are among the "best" overall measures of model performance [Willmott, 1982, Webster andOliver 2001].RMSE provides a measure of error size, but it is sensitive to outliers as it places a lot of weight on large errors [Chai and Draxler 2014].MAE provides an absolute measure of the size of the error and it is less sensitive to extreme values [Delbari & al. 2013]  Cross-validation of the model are shown in figures 6.a, 6.b and 6.c for the variables AB40, AB80 and AB120.Because the RMSE are closer to zero, fitted Exponential models yield better results and are adopted in this study.
After analysis of spatial dependence, we evaluated sites that had not been sampled in order to produce surface maps by ordinary kriging.(2) For AB=120 (figure 7.c) the disturbances as detected from surface measurements, which present depth higher than 20 m, seem to be not very significant regarding to surface study area.Thus, the majority of anomalies affect only the subsurface of the phosphate series

Figure 1 :
Figure 1: The basic principles of Schlumberger deviceProfiling resistivity surveys can be employed to detect lateral variations in resistivity.This type of survey provides estimates of the spatial variation in resistivity at some fixed electrode spacing[Chapellier 2000;Dabas and Caraire 2012].To determine the intrinsic resistivity for the different terms of the phosphatic series, we carried out, during the geophysical prospection in a parcel of 50 ha,5151 resistivity measurements as horizontal profiling (SYSCAL-R2 resistivity instrument IRIS Instruments) using the wellknown Schlumberger array, in order to map the spatial distribution of the sterile hardpan inclusions.The 5151 stations of the resistivity survey is a compilation of 51 profiles spaced at 20 m.For every profile, there were 101 stations separated by 5 m using the well-known Schlumberger array.

Figure 2 :
Figure 2: Typical semivariogram and its components.The technique of variogram (or semi-variogram) quantifies the spatial variability of a regionalized variable, and provides the input parameters for the spatial interpolation of kriging[Krige, 1951;Webster and Oliver, 2001].Let be r(x) and r(x+h) the values of electrical resistivity sampled at.
framework and methodology of workThe Ouled Abdoun basin is the largest phosphate basin in Morocco.It's located about 100 km in south-east of Casablanca.The phosphate deposits of Ouled Abdoun area belongs to the western Moroccan Meseta, commonly considered being stable.The local sedimentary deposits resulting from a large transgression occurred in mid-Cretaceous.It consists of :marly limestone and gypsum of Cenomanian, Turonian white limestones, Senonian marl and yellow marly limestones, phosphatic series dated from Maastrichtian to Ypresian and Lutetian calcareous Thersitean slab.The Neogene continental deposits cover locally the marine series.

Figure 4
Figure 4.a: Histogram of the resistivity data for AB=40m

Figure 4
Figure 4.b: Histogram of the resistivity data for AB=80m

Figure 4
Figure 4.c: Histogram of the resistivity data for AB=120m Nugget effect C 1 = Sill -Nugget effect A 1 = Range 1 Spatial Dependence Degree 2 Moderate

Figure 5
Figure 5.a: experimental semivariogram and fitted models for AB=40m . The definitions of MAE, MBE and RMSE are as follows: xi) and z*(xi) are the observed and predicted value respectively.To see how predicted values match closely the observed values, the RMSE can be used: the smaller the RMSE the better.

Figure 6
Figure 6.a: Cross validation graphs for AB40 Figures 7.a, 7.b and 7.c show isoresistivity Maps carried out at the nodes of the square grid of 1km x 1km with the Exponential model of variogram.

Table . 1
: Summary statistics of random samples

Table 2 :
Summary details of fitted Exponential models