Evaluation of reclamation success in an open-pit coal mine using integrated soil physical, chemical and biological quality indicators

Ecological Indicators 103 (2019) 182–193

Contents lists available at ScienceDirect

Ecological Indicators journal homepage: www.elsevier.com/locate/ecolind

Original Articles

Evaluation of reclamation success in an open-pit coal mine using integrated soil physical, chemical and biological quality indicators

T

⁎

Yamileth Domínguez-Haydara, , Elena Velásquezb, Jesús Carmonaa, Patrick Lavellec, Luis F. Chavezd,2, Juan J. Jiméneze,1 a

Universidad del Atlántico, Carrera 30 Número 8-49, Puerto Colombia, Colombia Universidad Nacional de Colombia, Carrera 32 Chapinero, vía Candelaria, Palmira, Colombia c IEES, Sorbonne University, Paris, France d Corporación Colombiana de Investigación Agropecuaria – Agrosavia, Centro de Investigación Palmira, Valle del Cauca, Colombia e Department of Biodiversity Conservation and Ecosystem Restoration, Pyrenean Institute of Ecology (IPE-CSIC), ES-22700 Jaca, Huesca, Spain b

A R T I C LE I N FO

A B S T R A C T

Keywords: Soil quality indices Dry tropical forest Soil macroinvertebrates Vegetation cover Technosols

The evaluation and success of ecological reclamation can be assessed by measuring soil physical, chemical and biological variables, either in isolation or combined into composite indicators. In this study we tested the suitability of biological, chemical and physical quality indicators–and their combination in a General Indicator of Soil Quality (GISQ)—to monitor soil quality in restored areas of an open cast coal mine. These composite indicators were computed with principal component analyses (PCA) and co-inertia multivariate analyses (CoIA). Our biological indicator showed a significant recovery of soil invertebrate communities along the chronosequence. Taxonomic richness increased from 7 taxa the 1st year to 13–17 taxa in 6–20-y sites, far less than in forests (19–20). Soil pH, bulk density and proportion of physical and macro aggregates were highest in the 1-y site, while soil organic matter (SOM), total N and proportion of biogenic aggregates were highest in older sites. In general, the three sub-indicators and the GISQ yielded the lowest values in the 1-y site (from 0.1 to 0.3 on a scale of 0.1 to 1.0), intermediate in the 16- and 20-y sites, (0.4–0.7) and reached highest values in the two forests (0.4–1.0). The GISQ methodology proved efficient in assessing progress in the reclamation process which should not be achieved only by monitoring vegetation cover changes.

1. Introduction By destroying the original plant cover and removing the soil, opencast mining has profound environmental impacts that need to be corrected once operations have been completed (Cooke and Johnson 2002; Nicolau and Moreno-de las Heras, 2005). Restoration of the natural capital in areas affected by mining activities is a necessary action to mitigate the effect of highly aggressive practices (Aronson et al., 2007). Restoration may take different forms–rehabilitation or reclamation–depending on the objectives set (SER, 2004; Lima et al., 2016). Assessment of ecosystem services in reclaimed mining sites is now widely performed, but the success of operations is still highly variable. Success is mainly evaluated based on the reconstitution of a dense and diverse plant cover (Ruiz-Jaen and Aide, 2005). However, several studies clearly show that soil and associated ecosystem services

restoration may still be very incipient in spite of a well-developed plant cover (Herrick et al., 2006). A healthy soil is necessary to sustain the development of vegetation and natural successional processes before a mature ecosystem is regenerated. Soil quality was defined by Doran and Parkin (1994) as “the capacity of soil to function, within ecosystem boundaries, to sustain biological productivity, maintain environmental quality, and promote plant and animal health”. In the context of mining operations, soil quality assessment may be complex because reclamation is performed on Technosols (IUSS Working Group WRB, 2014), i.e., a combination of mine spoils, natural topsoil and organic amendments with highly dynamic properties (Asensio et al., 2013; Capra et al., 2015). High- or low-quality soils are subjective expert judgments based on the kind of properties managers want soils to provide: high biodiversity in natural nutrient-poor soils of tropical forest, or chemical fertility and

⁎

Corresponding author. E-mail address: [email protected] (Y. Domínguez-Haydar). 1 ARAID Researcher. 2 Código 0000-0003-4251-9946. https://doi.org/10.1016/j.ecolind.2019.04.015 Received 10 November 2018; Received in revised form 7 March 2019; Accepted 5 April 2019 1470-160X/ © 2019 Published by Elsevier Ltd.

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

sampling period (October 2013–November 2014) a decrease in annual precipitation (407 mm) and an increase in evapotranspiration (2424 mm) was observed when compared with the historical trend of the last 10 years (Cerrejón, 2014). The Cerrejón coal mine is one of the largest in Latin America, with coal deposits within 69,000 ha of the territory. Operations in the mine were started in 1985 and an active reclamation program has been conducted in the last three decades. The reclamation of mined land begins even before any excavation is conducted. Topsoil is removed from forest areas where mining activities will start and further deposited in areas where active is implemented. A surface layer 30 cm deep of this natural topsoil was spread after a topographical reconstruction was performed. This substrate is denominated Technosol (IUSS Working Group WRB, 2015). At the mine, the official reclamation process is initiated by establishing one plot each year, following the procedure described above. The buffel grass (Cenchrus ciliaris L.) is sown in the plot as a technique to stabilize the soil (common procedure developed by technicians in the mine site). This was the starting point in our chronosequence of study sites. Two years later, 15–20 native species of the dry tropical forest, were planted in addition with the initial plant seedlings that derived from the soil seed bank of the natural topsoil (Domínguez-Haydar and Armbrecht, 2011; Gualdrón, 2011). No further actions were performed, except allowing regrowth of natural vegetation. Six rehabilitated areas were selected following a chronosequence of 1, 6, 8, 9, 16 and 20 years (-y). This approach is commonly used in studies conducted in rehabilitated mine areas where regular monitoring at plot level is not performed and true replicates rarely exist, for operational reasons (Domínguez-Haydar and Armbrecht, 2011; Mukhopadhyay et al., 2014; Ciarkowska et al., 2016). The dump material and soil cover used had similar origin and plants used were the same, creating a rather homogeneous condition at the onset of reclamation. Herbaceous plants from natural regeneration and OM debris, like branches and wood, were already observed in the 1-y site. In the 6y site, high densities of buffel grass and scattered trees (3-m tall) of Mimosa arenosa (Willd.) Poir. were present, and in the 8- and 9-y sites high densities of grass and trees with large canopies (4–5 m) were observed. At the 16 and 20-y site, the grass layer had decreased and leaf litter and bare soil covered a great proportion of the soil surface. Additionally, two patches of forests (Forest 1 and Forest 2) which resulted from natural regeneration (more than 40 years old), were selected as reference areas in the absence of primary forests that no longer exist in the region (Domínguez-Haydar and Armbrecht, 2011). Naturally regenerated areas with same ages as our reclamation sites did not either exist. For this reason, all comparisons were made among the different steps of the chronosequence and the two forests present at the site:

stable macro aggregation in cropped soils. Ecosystem restoration sets a different background for soil quality perception than, for example, agriculture or Technosol environments. In the extensive literature dedicated to soil quality indicators, three types of strategies have been used. One consists in simply measuring whatever biological, chemical or physical soil parameters are considered to be important indicators of soil properties. This may result in the need to measure long lists of parameters. Analysis of these parameters often shows that some parameters are variable and influenced by soil management options, while others are much more stable (Turbé et al., 2010; Askari and Holden, 2015). This observation has led to proposing limiting to a minimum set of relevant variables. Comprehensive studies of soil quality parameters show that many of them do covariate; as for example, pH with CEC, Ca, Mg and other bases, and with Al saturation or clay contents, in tropical forest or savanna soils (Grimaldi et al., 2014; Lavelle et al., 2014). Synthetic indicators that integrate these covariations have been proposed to overcome this inconvenience (Büneman et al., 2018). In this case, multivariate analyses and other numerical approaches allow identifying the main factors (= “latent variables”) that differentiate soils and define their quality (Andrews et al., 2004; Shukla et al., 2006; Janvier et al., 2007, Velásquez et al., 2007a; Turbé et al., 2010, De Paful Obade and Lal, 2016). In this study, the general indicator of soil quality (GISQ), a synthetic indicator of physical, chemical and biological quality was calculated with the methodology developed by Velásquez et al. (2007a,b). This index has been widely used to assess the effect of management options and general environmental changes–deforestation, and land use intensification on soil quality (Velásquez et al., 2012; Rousseau et al., 2012; Lavelle et al., 2014). At the open cast mine site of Cerrejón, in northern Colombia, reclamation has been implemented for 25 years and the percentage of cover attained in the reforestation plots is now close to that of the original forest (Domínguez-Haydar et al., 2018). However, we do not know how soil ecological functions and biodiversity have been recuperated. Under natural conditions, it may take 100 years for a forest to grow to full maturity, while only one centimeter of soil would be created over the same amount of time (Lavelle and Spain, 2005). Similar differences are expected to occur in an ecosystem under rehabilitation. If it were the case, specific actions might be necessary to accelerate the reclamation of soil quality and assist naturally slow processes. Reclamation actions in mining areas have special importance in biodiversity hot spots as is the case in Colombia (Arbeláez-Cortés, 2013). We evaluated soil quality indicators in a chronosequence (“space for time” approach, Frouz et al., 2013a) of rehabilitation plots aged 1–20 years. We tested the following hypotheses: (1) values of the physical, chemical and biological indicators and GISQ increase linearly with the age of reclamation while (2) the recovery of biological quality would happen sooner than that of physical and chemical properties in response to plant cover.

(a) Forest 1: Located close to the restored areas, this forest represented the initial pre disturbed condition and was the main source of the topsoil used for reclamation (“Caracolí” zone). In this forest, the subxerophytic vegetation was dominated by Prosopis juliflora (Sw.), Bulnesia arborea (Jacq) Engl., Platymiscium pinnatum (Jacq) Dugand, Cereus repandus (L.), Tabebuia billbergii (Bureau and K. Schum.), Bourreria exsucca (L). Jacq., Astronium graveolens Jacq., Guazuma ulmifolia Lam., Mimosa arenosa and the liana Serjania sp. (b) Forest 2 was a compensation area, excluded from mining activities. The tree layer was dominated by Aspidosperma polyneuron Müll.Arg. and Hura crepitans L. The shrub layer was represented by Cordia alba (Jacq.) Roem. and Schult., Machaerium sp. and Mimosa arenosa; some exotic African grasses like Andropogon sp., Brachiaria sp. and Panicum sp. were occasionally present (Gualdrón, 2011).

2. Materials and methods 2.1. Study area The study was carried out at the Cerrejón coal mine (11°8′56.34″N–10°57′36.16″N and 72°31′29.66″W–72°45′9.48″W), La Guajira department, Colombia, between 200 and 240 masl. (Fig. 1). The natural ecosystem in the area is a seasonally dry tropical forest (DTF), characterized by sub-xerophytic vegetation. The soil is classified as an Entisol (IGAC, 2014). Dry tropical forests globally represent 42% of all tropical forests (Brown and Lugo, 1982) and they are among the most endangered ecosystems in the world. Precipitation in the area is bimodal (March-April and October-November) with a yearly average value of 800 mm; mean annual temperature is 27.5 °C. During the

2.2. Experimental design In each area of the chronosequence and in the two control forests, four permanent monitoring plots (10 × 25 m) were selected. The 183

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Fig. 1. Location of the forest and rehabilitated areas included in this study.

different from classical methods that measure soil aggregate stability and size distribution. It was specifically designed to discriminate biogenic aggregates produced by soil organisms, roots and “physical” aggregates generated by physical processes where the biotic origin cannot be clearly identified. Soil texture was measured with the hydrometric method (Bouyoucos, 1962) and bulk density with the core method (Blake and Hartge, 1986). (c) Chemical sub-indicator – Soil chemical quality was assessed with eight selected variables. Total phosphorous (P) was measured by UV–VIS L-ascorbic spectrophotometry; Ca++, Mg++ and K+ by atomic absorption (ammonium acetate 1 N); pH (1:2 v/v water solution); organic matter (OM) content with the Walkley and Black volumetric method (Walkley and Black, 1934); total nitrogen with the Kjeldahl method; and effective cation exchange capacity (ECEC) as the sum of cations extracted with ammonium acetate 1 N (Burt, 2004).

minimum distance between plots was 200 m, and two sampling points were placed on opposite corners at 1-m distances outside each plot. In the two forests, plots were randomly selected following the same design. This resulted in eight sampling points per selected area along the chronosequence. Samples were taken at the end of the rainy season in October 2013 and October-November 2014. 2.3. Sampling protocol. Soil fauna survey and soil sampling. We evaluated soil quality/soil health along the chronosequence using three sub-indicators each one calculated from a different set of soil variables: (a) Biological sub-indicator – We used macroinvertebrate communities as an indicator of soil biodiversity and biological processes. As ecosystem engineers, invertebrates determine the composition and activity of microbial communities which justify their indicator status (Lavelle and Spain, 2005). A 25 × 25 × 20 cm3 soil core was taken in each sampling point of the plot and hand-sorted for macroinvertebrates (=all invertebrates visible to the naked eye), following standard ISO 23611-5:2011). Broad taxonomic units for soil macrofauna, either Order, Class or Subclass were used following identification keys in Borror et al. (1976). (b) Physical sub-indicator associated texture, bulk density and soil morphology analysis following Velásquez et al. (2007b) – Soil morphology: a 10 × 10 × 10 cm3 soil core was taken close to the soil macrofauna excavated pit, and carefully separated, following the procedure described in Velásquez et al. (2007a), into: (i) biogenic aggregates (> 5 mm), produced by macroinvertebrates recognized for their rounded forms, (ii) non-macroaggregated soil, (iii) root aggregates stuck to the roots (iv), physical aggregates (> 5 mm) produced by physical processes, with angular shapes and a smoother structure and, (v) organic debris, consisting of roots, leaf and woody debris and other materials of organic origin (organic debris). This protocol focused on the origin of aggregates

Additionally, plant cover was evaluated as a percent of canopy cover in each sampling point. This was measured by a Spherical Convex Crown Densiometer (Forestry Suppliers Inc., Jackson, USA) at 1 m from ground level (Lemmon, 1956). 2.4. Data analysis A four-step procedure was implemented to calculate the indicators in order to obtain each specific sub-indicator, i.e., biological, physical and chemical (ranging from 0.1 to 1) to be later combined in the General Indicator of Soil Quality (GISQ) following Velásquez et al. (2007a): (a) Exploration of data sets. A Principal Component Analysis (PCA) was performed on each data set, i.e., biological, physical and chemical. Coinertia analysis (CoIA) among pairs of data sets then allowed analyzing their covariations (Dolédec and Chessel, 1994; Dray et al., 2003). A matrix correlation coefficient (RV) was measured 184

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

(Table 1). Forest 2 had the richest community while the poorest community was found in the 1-y site, located at the opposite side of axis I. Axis II (10% of total variance explained) showed different taxa that characterized some specific sites. Lithobiomorpha (Chilopoda), Isoptera and Isopoda were associated along this axis with the 6-y, 8-y and 9-y sites, while Embioptera were specific to Forests and older sites. The lowest values of the biological sub-indicator were obtained in the 1-yr site. The indicator then increased and reached a stable value after 6–9 years, not different from the ones observed in the forests (Fig. 3a). In terms on composition (presence or absent of orders), the results of NMDS were stronger (stress 0.01) to demonstrate differences in the chronosequence. This showed a wide separation between areas of 1-y of rest, the intermediate sites conformed a group while that 16-y and 20-y remained grouped closer to forests (Fig. 4).

and its significance tested with a Monte Carlo test. Projection in a common factorial space of variables and samples from each respective data set enabled visualizing which specific variables covariate. We used the ade4 package in R. (b) Selection of variables to be used for GISQ computation. We selected the variables from each biological, chemical and physical data set that best differentiated the sites. Variables with the highest contribution to each factor eigenvalue (at least equal to half the maximum value calculated for variables on both axes) were selected. This allowed designing indicators with a reduced set of variables to be measured at further dates of the monitoring process. (c) Calculation of the different indicators. Each of the three sub-indicators of soil quality at a given site of the chronosequence was calculated following the procedure proposed by Velásquez et al. (2007a). (d) Lastly, a PCA was performed on a new data matrix that contained the values of the three sub-indicators and a GISQ was calculated following Velásquez et al. (2007a,b).

3.2. Soil physical quality variables and aggregates Soil texture was relatively homogeneous across sites, and no significant difference was observed in the value of silt and sand across the different sites. Average clay, sand and silt contents ranged from 21.5–28.5%, 39.9–46.5% and 32.9–36.1%, respectively. Soil bulk density was slightly higher in all the restored sites along the chronosequence compared with the reference forests (Table 2). A clear positive trend of an increase of the percentage of biogenic aggregates and concomitant decrease of the proportion of physical aggregates was observed. The first two axes of the PCA explained 39.7% of total variance (Fig. 5). The first axis of the PCA (22.2% variance explained) separated the less aged restored sites, where proportion of physical aggregates was highest. The second axis (17.5% of variability explained) separated the forests and the 20-y site (with high proportions of biogenic) from the most recently restored sites (1- and 6-y) which had the largest percentage of physical aggregates and low OM content (Table 2). As a result, the sub-indicator of physical quality exhibited lowest values in the 1-y and 6-y sites (Fig. 3b). The Montecarlo randomization test showed that the multivariate ordination was significant (p < 0.001, 25% explained variance).

The PCA and CoIA were performed with the ade4 library (Chessel et al., 2004), and the nlme library (Pinheiro et al., 2016) for linear mixed models computation in the R software (R Core Team, 2018). Non-metric multi-dimensional scaling (NMDS) based on the Jaccard index was used to examine the similarity between soil macrofauna community assemblages in the different sites. Canopy cover and chemical, physical and morphology variables were compared among areas with a non-parametric Kruskall–Wallis analysis of variance (ANOVA) and subsequently post-hoc Dunn test (Bonferroni adjustment at α < 0.05). Additionally, box-and-whisker plots were done for each sub-indicator across sites. Finally, we explored the relationship between GISQ and canopy cover with linear mixed models. In our study, the factor site (year after reclamation and forest type) was considered random to avoid violating the assumption of independence resulting from sets of correlated measurements (Chaves, 2010), and canopy cover was interpolated with natural cubic spline function with one degree of freedom (Welham et al., 2007; Wang, 2011). A model with random intercept and random slope was computed with lme4 library (Bates et al., 2015) for linear mixed models computation in the R software (R Core Team, 2018).

3.3. Chemical indicator variables 3. Results Soil pH ranged from 7.5 to 8.2, with relatively high nutrient contents, that increased with reclamation time. The forests showed the highest values of OM and N, and the lowest pH; significant differences were observed for these variables only between forests and the 1-y site (Table 2, Fig. 3c). With regards to the PCA (Fig. 6), the first axis explained 36.0% of the total variance, and sites were distributed along axis I depending on their chemical fertility, specifically for K+ contents and ECEC. Axis II explained 24.5% of the total variance and sites were chronologically separated, from 1- to 16- and 20-y and forests. The latter had the highest contents of OM and N while recently rehabilitated sites had high Ca++, Mg++ contents and ECEC, and the highest pH values (Table 2). The ordination was highly significant (Montecarlo randomization, p < 0.001), however the variance explained was relatively low (13%).

3.1. Biological indicator Soil macro invertebrates belonging to 18 broad taxonomic units (TU) were found. Forest 2 showed the highest taxa richness and density: 18 TU, 434 ind.m-2 (Table 1), while the lowest richness and density were observed in the 1-y site: 7 TU, 28 ind.m-2. In the rest of the chronosequence, TU ranged from 13 to 17, with the highest density observed in the 6-y site, i.e., 380 ind.m-2 and 8-y and 9-y sites (Table 1). Coleoptera, Hymenoptera and Isoptera presented the highest densities and were present along the chronosequence. Other groups with high densities were: Diplododa (Polydesmida), Isopoda and Zygentoma. The 1-y site was dominated by Coleoptera, with Carabidae as the most abundant family. The density of Arachnida and Hymenoptera, especially Formicidae, further increased in the 6-y and 8-y sites, while Embioptera was exclusive of the 16-y and 20-y sites and of both forests. Diplopoda (Polydesmida) were found in all sites between 9-y and 20-y and in Forest 2. Oligochaeta density was highest in Forest 1 and in the 8-y site (14 ind.m-2 and 5 ind.m-2) and were absent in the 1-y site and Forest 2. No group was exclusive of both forests (Table 1). Only two axes were retained from the PCA which explained 34% of the total variability (Fig. 2). The permutation test showed significant differences among sites (p < 0.01). Axis I (24% of total variability) separated sites according to the number of taxa and macrofauna density

3.4. Changes in the General Soil Quality Index In general, all sub-indicators and the GISQ noticeably increased with the age of the restored site (Fig. 3d) Average values of the soil biological macrofauna sub-indicator ranged from a minimum of 0.19 at the 1-y site to 0.36 and 0.47 in forest 1 and 2, respectively (Table 3). The 6-y to 20-y restored sites had intermediate values, i.e. from 0.27 to 0.33, and did not exhibit a clear temporal pattern. The highest recovery was observed in the physical sub-indicator (Fig. 3d). The CoIA between soil macrofauna and chemical variables showed 185

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Table 1 Mean densities (ind m−2) and standard error (italic) of soil macrofauna along the chronosequence of rehabilitated areas and in the undisturbed forests. 1y

6y

8y

9y

16 y

20 y

Forest 1

Forest 2

Acari

Mean SE Range

0.00 0.00 0.00

0.00 0.00 0.00

0.00 0.00 0.00

0.00 0.00 0.00

1.00 1.00 0–16

0.00 0.00 0.00

4.00 1.79 0–16

2.00 2.00 0–32

Aranae

Mean SE Range

2.00 5.47 0–16

4.00 9.24 0–32

9.00 17.50 0–48

0.00 0.00 0.00

3.00 6.45 0–16

3.00 8.70 0–32

6.00 9.91 0–32

16.00 18.48 0–48

Coleoptera

Mean SE Range

13.00 4.67 0–64

40.00 13.15 0–144

22.00 6.83 0–96

12.00 5.16 0–64

16.00 7.30 0–112

23.00 6.36 0–64

33.00 11.73 176.00

24.00 8.00 112.00

Coleoptera Larvae

Mean SE Range

3.00 2.18 0–32

5.00 2.41 0–32

12.00 3.10 0–32

30.00 22.34 0–352

13.00 3.34 0–32

6.00 3.54 0–48

17.00 6.10 64.00

37.00 9.09 96.00

Embioptera

Mean SE Range

0.00 0.00 0–0

0.00 0.00 0–0

0.00 0.00 0–0

0.00 0.00 0–0

7.00 4.12 0–48

7.00 2.52 0–32

11.00 5.60 80.00

14.00 7.57 96.00

Escolopendro-morpha

Mean SE Range

0.00 0.00 0.00

0.00 0.00 0.00

3.00 2.18 0–32

2.00 1.37 0–16

5.00 1.91 0–16

2.00 1.37 0–16

3.00 1.61 0–16

5.00 2.41 0–32

Gastropoda

Mean SE Range

0.00 16.00 0–16

4.00 1.79 0–16

7.00 3.57 0–48

10.00 3.54 0–48

8.00 2.92 0–32

4.00 2.31 0–32

3.00 2.18 0–32

9.00 3.26 0–32

Hymenoptera

Mean SE Range

3.00 2.18 0–32

247.00 134.43 0–1952

100.00 75.32 0–1216

50.00 24.65 0–336

40.00 18.70 0–304

11.00 6.96 0–112

98.00 77.01 1232.00

166.00 105.40 1696.00

Isoptera

Mean SE Range

1.00 1.00 0–16

37.00 21.94 0–320

17.00 10.48 0–160

61.00 37.87 0–464

22.00 15.10 0–192

22.00 10.42 0–144

20.00 14.35 208.00

70.00 28.65 352.00

Isopoda

Mean SE Range

0.00 0.00 0–0

18.00 6.67 0–80

48.00 24.92 0–352

11.00 5.97 0–96

3.00 2.18 0–32

2.00 1.37 0–16

11.00 11.00 176.00

34.00 22.62 368.00

Julida

Mean SE Range

0.00 0.00 0

0.00 0.00 0

0.00 0.00 0

5.00 5.00 0–80

0.00 0.00 0.00

0.00 0.00 0.00

5.00 5.00 0–80

2.00 1.37 0–16

Lepidoptera larvae

Mean SE Range

4.00 1.79 0–16

4.00 2.31 0–32

2.00 1.37 0–16

8.00 3.86 0–48

1.00 1.00 0–16

3.00 2.18 0–32

2.00 1.37 0.00

4.00 2.31 0–32

Geophilomorpha

Mean SE Range

2.00 1.37 0–16

7.00 2.91 0–32

3.00 2.18 0–32

2.00 1.37 0–16

0.00 0.00 0.00

1.00 1.00 16.00

0.00 0.00 0.00

0.00 0.00 0.00

Neuroptera

Mean SE Range

0.00 0.00 0.00

1.00 1.00 0–16

1.00 1.00 0–16

1.00 1.00 0–16

0.00 0.00 0.00

2.00 1.37 0–16

2.00 1.37 0–16

0.00 0.00 0.00

Orthoptera

Mean SE Range

0.00 0.00 0.00

2.00 1.37 0–16

1.00 1.00 0–16

1.00 1.00 0–16

0.00 0.00 0.00

0.00 0.00 0.00

0.00 0.00 0.00

5.00 2.41 0–32

Oligochaeta

Mean SE Range

0.00 0.00 0–0

3.00 2.18 0–32

5.00 3.17 0–48

3.00 2.18 0–32

3.00 2.18 0–32

4.00 3.10 0–48

0.00 0.00 0.00

14.00 4.35 48.00

Polidesmida

Mean SE Range

0.00 0.00 0.00

0.00 0.00 0.00

0.00 0.00 0.00

2.00 1.37 0–16

6.00 2.00 0–16

11.00 2.82 0–32

1.00 1.00 0–16

7.00 2.52 0–32

Pseudoescorpio-nidae

Mean SE Range

0.00 0.00 0.00

0.00 0.00 0.00

1.00 1.00 0–16

1.00 1.00 0–16

0.00 0.00 0.00

0.00 0.00 0.00

5.00 1.91 0–16

6.00 2.88 0–32

Zygentoma

Mean SE Range

0.00 0.00 0–0

12.00 6.11 0–80

38.00 10.62 0–144

21.00 10.71 0–128

12.00 4.73 0–64

9.00 6.02 0–96

4.00 2.31 0.00

8.00 3.27 0–48

7.00 26.00

13.00 380.00

15.00 260.00

17.00 228.00

15.00 148.00

17.00 144.00

14.00 209.00

18.00 423.00

Richness§ Density§§

SE = standard error. § Number of orders per site. §§ Number of individuals m−2.

186

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Fig. 2. Principal component analysis (PCA) biplot of macroinvertebrate communities with the orders and the projection of forests and rehabilitated areas as objects in the plane defined by the first two factors. Rich = Number of orders encountered, L = Larvae, Arac = Arachnida, Cole = Coleoptera, L. Col = Coleoptera larvae, Embi = Embioptera, Esco = Escolopendromorpha, Gast = Gastropoda, Hyme = Hymenoptera, Isop = Isopoda, Isopt = Isoptera, L.Lep = Lepidoptera larvae, Lith = Lithibiomorpha, Hapl = Haplotaxida, Poly = Polydesmida, Zyg = Zygentoma.

Fig. 3. Box and whisker plots of sub-indicators and plant cover in each site evaluated.

The CoIA between soil physical and chemical variables (RV = 0.65, p < 0.05), showed a strong correlation between OM and N contents, and proportion of litter and biogenic aggregates (Fig. 7c). The first two axes of the PCA analysis of the three sub-indicators explained 81.4% of total variance (Fig. 8). All indicators projected on the same side of the first axis (F1: 55.2%). A logical temporal sequence of recovery was observed along the first axis. The effect of time illustrated by contrasted projection of the 1-yr site on the right-hand side, and forests and old sites on the left-hand side was highly significant with a rapid increase in the first 6 years (p < 0.001). The second axis (explained 26.2% of total variability) and clearly opposed the soil macrofauna sub-indicator to the chemical and physical ones. The canopy cover steadily increased with the age of the chronosequence, reaching a 50% coverage 6 years after the onset of the reclamation process (Fig. 3e). The 16- and 20-y sites had the highest canopy cover, similar to the forests (Table 2). A positive relationship was observed between canopy cover and the GISQ (Fig. 9, Table 4). However, the relationship of canopy cover with GISQ varied within

strong but not significant covariation (RV = 0.60; p > 0.05). Axis I explained 78.2% of total variability and separated both forests (Fig. 7a) from sites along the chronosequence, showing a positive relationship between total abundance and richness of communities, and OM and N contents. The lowest abundance and richness of macrofauna was obtained at the early stages of reclamation, where soil pH, Ca++ and Mg++ contents were highest. Embioptera were only present in forests and in 16- and 20-y sites. Axis II (16% of total variance explained) separated forests from the rest of sites along the chronosequence. Soil macrofauna communities also showed a high covariation with soil physical properties (RV = 0.60; p < 0.05), with an ordination roughly following the chronosequence (Fig. 7b). The first axis (59% of the total variance) separated the forests which had a higher proportion of biogenic aggregates and OM content–clearly associated to the highest density of macrofauna groups. Root aggregates were strongly related to the 8-y site. Total richness and abundance of soil macrofauna was lowest in the 1-y site, which was characterized by a larger proportion of physical aggregates and higher bulk density than the forests. 187

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Fig. 4. Non-metric multidimensional scaling (NMDS) Multidimensional scaling output of ant guilds in sites with different rehabilitation ages and forests (F) at Cerrejón coal mine, Guajira, Colombia.

Table 2 Mean values and standard error (italic) of chemical, morphology and physical soil variables and plant cover in rehabilitated areas and in the undisturbed forests. Kruskall Wallis test was performed and post hoc comparisons using the Dunn test. Letters associated with mean values indicate significant differences between means of P < 0.05; the absence of letters indicates that no significant differences were found. Variables/Age (y) Chemical variables pH P (ppm) Ca (mq) Mg (meq) K (meq) ECEC OM (%) N (%) Physical variables SNA (gr) BA (gr) PA (gr) RA (gr) OD (gr) St (gr) BD (gr) S (%) Si (%) C (%) Canopy cover CC (%)

1

6

8

9

16

20

Forest-1

Forest-2

8.18 a 0.15 44.86 ab 25.73 55.62 15.1 6.2 1.89 0.94 0.53 59.05 21.76 2.24 a 0.91 0.11 a 0.05

8.15 ab 0.31 42.73 ab 26.41 53.95 13.65 5.33 2.42 1.29 0.83 61.3 14.25 2.71 abc 0.88 0.13 abc 0.04

8.08 ac 0.22 40.36 ab 25.23 50.12 12.11 5.28 2.42 1.09 0.82 54.09 18.15 2.58 a 0.99 0.13 ac 0.04

8.17 a 0.2 52.82 b 35.65 55.27 19.16 5.21 2.19 1.51 1.27 51.24 20.05 2.77 a 2.7 0.14 a 0.13

8.17 ab 0.16 20.32 a 21.13 54.68 18.15 4.85 1.7 1.44 1.24 62.61 17.12 3.06 abc 1.18 0.16 abc 0.05

8.13 abc 0.23 31.31 ab 21.65 51.7 13.98 5.34 3.39 1.3 1.11 58.28 17.03 3.05 abc 1.1 0.15 abc 0.06

7.78 b 0.58 17.68 a 11.42 44.53 20.76 3.94 1.19 1.12 1.18 59 16.61 4.08 b 1.24 0.21 b 0.05

7.46 c 0.87 20.89 ab 17.33 45.2 15.91 4.75 1.97 1.61 2.06 58.82 18.38 3.93 c 1.12 0.2 c 0.05

491.63 191.64 0a 0 829.91 a 211.03 1.12 ac 3.26 1.23 a 1.62 15.7 ab 14.4 1.25 a 0.1 40.37 5.71 34 6.06 25.62 ab 4.08

564.46 172.85 312.78 abc 239.12 196.94 b 258.66 1.72 ac 3.09 2.53 ab 3.17 17.33 ab 24.57 1.1 ab 0.17 39.87 7.35 31.87 5.58 28.25 b 5.15

590.4 162.78 586.69 bcd 154.8 10.99 bc 17.21 2.83 ac 6.78 3.06 a 4.56 18.32 ab 26.15 1.26 a 0.12 45.87 5.23 30 4.67 24.12 ab 4.64

557.82 197.18 628.42 bcd 283.21 22.31 bc 42.47 2.11 ac 3.25 2.31 a 2 17.62 ab 22.98 1.19 ab 0.09 43.62 5.01 33.5 5.48 22.87 ab 3.42

572.71 133.59 608.8 bcd 194.37 10.57 bc 14.86 0.39 ac 1.08 2.08 ab 1.52 64.46 b 87.67 1.17 ab 0.1 46.12 5.53 30.25 6.01 23.62 ab 2.65

446.4 140.02 696.53 d 180.93 21.13 bc 55.92 0c 0 1.8 ab 1.77 24.63 ab 41.58 1.2 ac 0.06 42.5 4.16 33 4.84 24.5 ab 3.54

534 196.63 515.53 cd 225.81 0c 0 4.1 b 6.66 4.14 b 3.4 13.15 a 30.52 1.09 b 0.09 46.5 10.64 32 6.85 21.5 ab 7.28

477.86 195.16 586.65 cd 176.84 0.04 c 0.17 1.91 a 2.66 4.18 b 2.76 14.7 a 36.38 1.05 bc 0.09 41.12 10.6 36.12 6.86 22.75 ab 6.56

0.16 a 0

52.81 ab 27.1

55.54 ab 31.56

75.3 bc 31.01

94.15 d 10.32

93.11 dc 9.26

92.2 bd 8.15

99.61 d 0.77

ECEC = Effective cation exchange capacity, OM = Organic matter, NAS = Non aggregate soil, BA = Biogenic aggregate, PA = Physical aggregate, RA = Root aggregate, OD = Organic debris, St = Stones, BD = Bulk density, S = Sand, Si = Silt, C = Clay, CC = Canopy cover. 188

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Fig. 5. Principal component analysis (PCA) biplot of physical and morphological variables and projection of forests and rehabilitated areas in the plane defined by factors 1 and 2. NAS = Non soil aggregate, BA = Biogenic aggregate, PA = Physical aggregate, RA = Root aggregate, OD = Organic debris, St = Stones, BD = Bulk density, S = Sand, Si = Silt, C = Clay.

sites. Recent 6- to 8-y old sites with 100% canopy cover had lower values for the GISQ than the older sites (Fig. 9).

Table 3 Mean values of sub-indicators and General Soil Quality Index (GISQ) of 16 points per site (eight points for each sample). Sub-indicators

4. Discussion Soil quality increased along the chronosequence, although differently depending on the sub indicator considered. PCAs performed on the three different data matrices provided similar separation of chronosequence sites, roughly ranking these based on their ages – from 1- to 20-y-old reclamation plots and forests. The GISQ indicator integrates sub-indicators that were calculated from the corresponding soil quality variables that are linked to soil biological, chemical, and physical processes involved (Velásquez et al., 2007a). In this study we aimed to test the use of the GISQ methodology in conditions where a relatively small set of soil variables was used when compared with previous studies (Velásquez et al., 2007a; Grimaldi et al., 2014). This general index of soil quality increased with time of reclamation; however, the dynamics of each of the sub indicators showed that their progress did not keep the same temporal pace. The physical indicator reached its highest value after 8 years, while the biological indicator between 8 and 16 years, and the chemical indicator did not change.

GISQ

Age/Forest

Biologicalψ

Chemical

Physical

Mean (range)

1-y 6-y 8-y 9-y 16-y 20-y Forest 1 Forest 2

0.19 0.30 0.33 0.27 0.28 0.27 0.36 0.47

0.26 0.31 0.31 0.30 0.33 0.33 0.43 0.46

0.19 0.55 0.76 0.76 0.76 0.79 0.76 0.78

0.18 0.42 0.52 0.53 0.50 0.50 0.55 0.62

ψ

(0.1–0.3) (0.2–0.5) (0.4–0.6) (0.3–0.7) (0.4–0.6) (0.4–0.7) (0.5–0.7) (0.4–1.0)

Soil macrofauna groups.

age. The soil macrofauna sub-indicator increased along the chronosequence, to a maximum of 0.36–0.47 in forests. Soil macrofauna density recovered quickly with high values at year 6, while maximum richness was recorded at year 9. The recovery of soil organisms during reclamation depends on the gradual colonization from external species pools, release from disturbance and the response of the component species to changes in environmental conditions (Kardol et al., 2009). Local constraints may also add their effect on the composition of communities (Belyea and Lancaster, 1999). Although not measured in this study, an important factor affecting colonization is species dispersal ability. The vicinity of the selected areas of the chronosequence to

4.1. Soil macrofauna composition along the chronosequence Soil macroinvertebrate communities showed a progressive recovery pattern, with the largest part of the biodiversity recovered after 6 years. NMDS analysis showed that macrofauna responded to rehabilitation

Fig. 6. Principal component analysis (PCA) biplot of chemical variables and projection of forests and rehabilitated areas in the plane defined by factors 1 and 2. ECEC = Effective cation exchange capacity, OM = Organic matter. 189

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Fig. 7. Analysis of co-variations using the co-inertia analysis (CoIA). (a) Between soil macrofauna and chemical variables. P > 0.05, RV = 0.6; (b) Between soil macrofauna and physical and morphological variables. P < 0.05, RV = 0.6; and (c) Between chemical variables, and physical and morphological variables P < 0.05, RV = 0.65. The vectorial correlation coefficient (RV) is equivalent to the r-correlation coefficient (ranges between 0 and 1), and is a measure of similarity between squared symmetric matrices (Escoufier, 1973). Arac = Arachnida, Cole = Coleoptera, L.Col = Larvae of Coleoptera, Embi = Embioptera, Gast = Gastropoda, Hyme = Hymenoptera, Iso = Isopoda, Isop = Isoptera, L.Lepi = Lepidoptera (Larvae), Hapl = Haplotaxida, Zyn = Zyngentoma.

Fig. 8. Projection of forests and rehabilitated areas in factorial space defined by the Principal component analysis (PCA) of subindicators.

Barrera, 2007). In the 6-, 8- and 9-y sites, the presence of a herbaceous layer (especially C. ciliaris) and of trees created a more diverse microenvironment favoring high densities of macrofauna. Omnivore soil macrofauna like ants, specially the genera Pheidole and Solenopsis, and families Cucujidae and Staphylinidae (Coleoptera) were abundant. On the contrary, in the 16- and 20-y sites, large areas of bare soil appeared

natural forests (1–2 km apart) may have eased the colonization of those areas through vegetation corridors (Domínguez-Haydar and Armbrecht, 2011; Frouz et al., 2013a). For example, at the 1-y site, ground beetles were the most abundant group; this could be favored by the proximity to Forest 1, a habitat source (Wanner and Dunger, 2002), and by the presence of other organisms feeding on plant and wood debris contained in the topsoil on which these Carabidae prey (Granados and 190

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

the 1-y site and the forest, the OM contents found in the former were higher in comparison with other mine sites from tropical areas (Mukhopadhyay et al., 2014). This legacy was demonstrated with soil micromorphology analysis of thin sections of resin-impregnated soil cores (Domínguez-Haydar et al., 2018). 4.3. Relation between GISQ and vegetation cover The indicator was positive and significantly correlated with vegetation cover along the chronosequence (Table 4). Canopy cover per se does not stand for a good indicator of recovery of soil quality, given that low values of GISQ were found in areas with high canopy cover. However, a good soil quality can be guaranteed in habitats with more complex vegetation structure. The response of plant-driven above- and below-ground processes have indeed a strong influence on soil quality (Kardol and Wardle, 2010; Menta et al., 2014; Mukhopadhyay et al., 2016). Fig. 9. Relation between mean values of the General Index of Soil Quality (GISQ) and plant cover.

5. Conclusion In our study the aim was to test the suitability of the GISQ methodology to assess minesite restoration success. Soil quality/soil health changes along the chronosequence were evidenced when all sub-indicators were integrated in the GISQ. One of the drawbacks of this study was the number of selected soil variables for the analysis, which was linked to the characteristics of the site and the possibility of performing soil analyses in the company’s premises. Despite this, the use of this specific soil quality indicator in the site-specific conditions has resulted in an efficient tool for a small number of variables in the data set. However, a major taxonomic resolution in highly diverse macrofauna groups like Coleoptera should be sought. By so doing, a more precise assessment of the soil macrofauna subindicator in tropical habitats would be achieved. Another important issue relates to the nature of soil, i.e., the mix of horizons. From our results we suggest to analyze in forthcoming studies the soil at different depths (0–5 cm, 5–10 cm, etc) in order to detect more discrete changes in the topsoil that is being “constructed”. In addition to other soil quality indicators, the advantage of GISQ, is its highly suitability and performance for soil quality assessment procedure in Technosols.

Table 4 Results of the linear mixed model when using forest cover as predictor variable for GISQ.

(Intercept) Cover §

Estimate§

Std.Error

t-value

p-value

0.266 0.364

0.066 0.103

4.022 3.524

0.0146 0.0149

Degrees of freedom = 54.

as a consequence of the reduction of grasses, which resulted in low soil macrofauna abundance. Only other specialized soil macrofauna groups, such as Embioptera and Polydesmida, i.e., litter detritivores related to sites with dense plant cover (Doblas-Miranda et al., 2008; Frouz et al., 2013b; Bogyó et al., 2015), and Pseudoscorpionida (predators) were present in the more aged restored areas. These may be contributing to the position of the 20-yr and 16-yr with respect to forests in the NMDS. However, soil macrofauna richness and composition at advanced stages along the chronosequence did not reach the values obtained in the forests.

Acknowledgements

4.2. Change of soil physical and chemical properties along the restored chronosequence

This study was supported by a COLCIENCIAS grant, Colombia (Code: 1116-569-34827) and “Vicerrectoría de Investigaciones, Extensión y Proyección Social” of Universidad del Atlántico. Particular thanks are extended to IPE-CSIC (Spain), to the “Departamento de Gestión Ambiental” of “Cerrejón” for logistic assistance, Y. Barros for field assistance, and K. Álvarez and L. Quijano for lab assistance. The ARAID, Spain foundation is acknowledged for support to J.J. Jiménez. We thank two anonymous reviewerstheir useful comments in a previous version of this paper.

The sub-indicator of soil physical quality increased quicker than the other two sub-indicators with the age of the restored site, reaching values close to the forest in a relatively short period. The correlation found in the CoIA between biogenic aggregates, OM and N contents with litter and soil macrofauna illustrated the bioturbation effect of the macrofauna as ecosystem engineers (Brussaard et al., 2007; Velásquez et al., 2012; Lavelle et al., 2016), which are attributes characteristic of forested sites. For example, although Oligochaeta abundance was low, their activity can be indirectly measured by the formation of biogenic aggregates (Velásquez et al., 2007b). This technique has been successfully used in other studies and later validated with near infrared spectroscopy (NIRS) that confirmed the different origins of aggregate classes based on their spectral signatures (Hedde et al., 2005; Velásquez et al., 2007a; Cunha et al., 2016; Zangerlé et al., 2016). The increase in the proportion of biogenic aggregates and the decrease of physical aggregates along the chronosequence demonstrated that soil structure was improved with the age of reclamation (Domínguez-Haydar et al., 2018). This stimulates pedogenic processes and evolution of these Technosols (Wilkinson et al., 2009; Leguédois et al., 2016). The topsoil has a chemistry legacy translocated from forests that sets the conditions for plant establishment and support of other soil functions. For example, despite differences observed in OM content between

References Andrews, S.S., Karlen, D.L., Cambardella, C.A., 2004. The soil management assessment framework. Soil Sci. Soc. Am. J. 68 (6), 1945–1962. https://doi.org/10.2136/ sssaj2004.1945. Arbeláez-Cortés, E., 2013. Knowledge of Colombian biodiversity: published and indexed. Biodivers. Conserv. 22, 2875–2906. https://doi.org/10.1007/s10531-013-0560-y. Aronson, J., Milton, S.J., Blignaut, J.N., 2007. Restoring Natural Capital: Science, Business, and Practice. Society for Ecological Restoration International. Island Press https://doi.org/10.1111/j.1526-100X.2008.00392.x. Asensio, V., Vega, F.A., Andrade, M.L., Covelo, E.F., 2013. Technosols made of wastes to improve physico-chemical characteristics of a copper mine soil. Pedosphere 23, 1–9. https://doi.org/10.1016/S1002-0160(12)60074-5. Askari, M.S., Holden, N.M., 2015. Quantitative soil quality indexing of temperate arable management systems. Soil Tillage Res. 150, 57–67. https://doi.org/10.1016/j.still. 2015.01.010Get.

191

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Oliveira, M.N., Ramirez, B., Rodríguez, G., Sarrazin, M., Da Silva, M.L., Costa, L.G.S., De Souza, S.L., Veiga, I., Velásquez, E., Lavelle, P., 2014. Ecosystem services of regulation and support in Amazonian pioneer fronts: searching for landscape drivers. Landscape Ecol. 29, 311–328. https://doi.org/10.1007/s10980-013-9981-y. Gualdrón, R., 2011. Hacia la rehabilitación de las tierras intervenidas por la minería a cielo abierto (accessed 02 January 2018) [in Spanish]. http://www.cerrejon.com/ wp-content/uploads/2017/12/Hacia-la-rehabiltacion-de-tierras-min.pdf. Hedde, M., Lavelle, P., Joffre, R., Jiménez, J.J., Decaëns, T., Jimenez, J.J., Decaens, T., 2005. Specific functional signature in soil macro-invertebrate biostructures. Funct. Ecol. 19 (5), 785–793 5 10.1111/j.1365-2435.2005.01026.x. Herrick, J.E., Schuman, G.E., Rango, A., 2006. Monitoring ecological processes for restoration projects. J. Nat. Conserv. 14, 161–171. https://doi.org/10.1016/j.jnc.2006. 05.001. Instituto Geográfico Agustín Codazzi (IGAC). 2014. Manejo de Suelos Colombianos. Bogotá. 323 p. v. IUSS Working Group WRB, 2014. World reference base for soil resources 2014. International soil classification system for naming soils and creating legends for soil maps, World Soil Resources Reports No. 106. FAO, Rome. IUSS Working Group WRB, 2015. World reference base for soil resources 2014, update 2015. In: International Soil Classification System for Naming Soils and Creating Legends for Soil Maps. World Soil Resources Reports No. 106. FAO, Rome. https:// doi.org/10.1017/S0014479706394902. Janvier, C., Villeneuve, F., Alabouvette, C., Edel-Hermann, V., Mateille, T., Steinberg, C., 2007. Soil health through soil disease suppression: which strategy from descriptors to indicators? Soil Biol. Biochem. 39, 1–23. https://doi.org/10.1016/j.soilbio.2006.07. 001. Kardol, P., Wardle, D.A., 2010. How understanding aboveground–belowground linkages can assist restoration ecology. Trends Ecol. Evol. 25, 670–679. https://doi.org/10. 1016/j.tree.2010.09.001. Kardol, P., Newton, J.S., Martijn Bezemer, T., Maraun, M., van der Putten, W., 2009. Contrasting diversity patterns for soil mites and nematodes in secondary succession. Acta Oecologica 35, 603–609. https://doi.org/10.1016/j.actao.2009.05.006. Lavelle, P., Spain, A., 2005. Soil Ecology, 2nd Ed. Kluwer Academic Publishers, Dordrecht, pp. 654. Lavelle, P., Rodríguez, N., Arguello, O., Bernal, J., Botero, C., Chaparro, P., Gómez, Y., Gutiérrez, A., Hurtado, M.P., Loaiza, S., Pullido, S.X., Rodríguez, E., Sanabria, C., Velásquez, E., Fonte, S.J., 2014. Soil ecosystem services and land use in the rapidly changing Orinoco River Basin of Colombia. Agric. Ecosyst. Environ. 185, 10. https:// doi.org/10.1016/j.agee.2013.12.0206-117. Lavelle, P., Spain, A., Blouin, M., Brown, G., Decaëns, T., Grimaldi, M., Jiménez, J.J., McKey, D., Mathieu, J., Velasquez, E., Zangerlé, A., 2016. Ecosystem engineers in a self-organized soil. Soil Sci. 181, 91–109. https://doi.org/10.1097/SS. 0000000000000155. Leguédois, S., Séré, G., Auclerc, A., Cortet, J., Huot, H., Ouvrard, S., Watteau, F., Schwartz, C., Morel, J.L., 2016. Modelling pedogenesis of Technosols. Geoderma 262, 199–212. https://doi.org/10.1016/j.geoderma.2015.08.008. Lemmon, P.E., 1956. A spherical densiometer for estimating forest overstory density. Forest Sci. 2, 314–320. Lima, A.T., Mitchell, K., O’Connell, D.W., Verhoeven, J., van Cappellen, P., 2016. The legacy of surface mining: remediation, restoration, reclamation and rehabilitation. Environ. Sci. Policy 66, 227–233. https://doi.org/10.1016/j.envsci.2016.07.011. Menta, C., Conti, F.D., Pinto, S., Leoni, A., Lozano-Fondón, C., 2014. Monitoring soil restoration in an open-pit mine in northern Italy. Appl. Soil Ecol. 83, 22–29. https:// doi.org/10.1016/j.apsoil.2013.07.013. Mukhopadhyay, S., Maiti, S.K., Masto, R.E., 2014. Development of mine soil quality index (MSQI) for evaluation of restoration success: a chronosequence study. Ecol. Eng. 71, 10–20. https://doi.org/10.1016/j.ecoleng.2014.07.001. Mukhopadhyay, S., Masto, R.E., Yadav, A., George, J., Ram, L.C., Shukla, S.P., 2016. Soil quality index for evaluation of reclaimed coal mine spoil. Sci. Total Environ. 542, 540–550. https://doi.org/10.1016/j.scitotenv.2015.10.035. Nicolau, J.M., Moreno-de las Heras, M., 2005. Opencast mining restoration. In: Mansourian, S., Vallauri, D., Dudley, N. (Eds.), Forest Restoration in Landscapes: Beyond Planting Trees. Springer, NY, pp. 370–378. https://doi.org/10.1007/0-38729112-1. Pinheiro, J., Bates, D., DebRoy, S., Sarkar, D., 2016. Linear and Nonlinear Mixed Effects Models [R package nlme version 3.1-128]. R Core Team, 2018. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria http://www.R-project.org/. Rousseau, G.X., Deheuvels, O., Rodriguez-Arias, I., Somarriba, E., 2012. Indicating soil quality in cacao-based agroforestry systems and old-growth forests: The potential of soil macrofauna assemblage. Ecol. Ind. 23, 535–543. https://doi.org/10.1016/j. ecolind.2012.05.008. Ruiz-Jaen, M.C., Aide, T.M., 2005. Restoration success: how is it being measured? Restor. Ecol. 13, 569–577. https://doi.org/10.1111/j.1526-100X.2005.00072.x. SER, 2004. The SER International primer on ecological Restoration. Version 2. Society for Ecological Restoration International Science and Policy Working Group. www.ser. org/pdf/primer3.pdf. Shukla, M.K., Lal, R., Ebinger, M., 2006. Determining soil quality indicators by factor analysis. Soil Tillage Res. 87, 194–204. https://doi.org/10.1016/j.still.2005.03.011. Turbé, A., De Toni, A., Benito, P., Lavelle, P., Lavelle, P., Ruiz, N., Mudgal, S., 2010. Soil biodiversity: functions, threats and tools for policy makers. Report for European Commission (DG Environment). Bio Intelligence service, IRD, NIOO https://halbioemco.ccsd.cnrs.fr/bioemco-00560420/. Velásquez, E., Fonte, S.J., Barot, S., Grimaldi, M., Desjardins, T., Lavelle, P., 2012. Soil macrofauna-mediated impacts of plant species composition on soil functioning in Amazonian pastures. Appl. Soil Ecol. 56, 43–50. https://doi.org/10.1016/j.apsoil.

Bates, D., Maechler, M., Bolker, B., Walker, S., 2015. Fitting linear mixed-effects models using lme4. J. Stat. Softw. 67 (1), 1–48. https://doi.org/10.18637/jss.v067.i01. Belyea, L.R., Lancaster, J., 1999. Assembly rules within a contingent ecology. Oikos 86, 402–416. https://doi.org/10.2307/3546646. Blake, G.R., Hartge, K.H., 1986. Bulk density. In: Klute, A. (Ed.), Methods of Soil Analysis. American Society of Agronomy, Madison, WI, pp. 363–375. Bogyó, D., Magura, T., Nagy, D.D., Tóthmérész, B., 2015. Distribution of millipedes (Myriapoda, Diplopoda) along a forest interior – forest edge – grassland habitat complex. Zookeys 510, 181–195. https://doi.org/10.3897/zookeys.510.8657. Borror, D., Delong, D., Triplehorn, C., 1976. An Introduction to the Study of Insects. Holt, Rinehart and Winston, New York, pp. 852. Bouyoucos, G.J., 1962. Hydrometer method improved for making particle size analyses of soils. Agron. J. 54, 464. https://doi.org/10.2134/agronj1962. 00021962005400050028x. Brown, S., Lugo, A.E., 1982. Storage and production of organic matter in tropical forests and their role in the global carbon cycle. Biotropica 14, 161–187. https://doi.org/10. 2307/2388024. Brussaard, L., Pulleman, M.M., Ouédraogo, É., Mando, A., Six, J., 2007. Soil fauna and soil function in the fabric of the food web. Pedobiologia 50, 447–462. https://doi. org/10.1016/j.pedobi.2006.10.007. Get rights and content. Bünemann, E.K., Bongiorno, G., Bai, Z., Creamer, R.E., De Deyn, G., de Goede, R., Fleskens, L., Geissen, V., Kuyper, T.W., Mäder, P., Pulleman, M., Sukkel, W., van Groenigen, J.W., Brussaard, L., 2018. Soil quality – a critical review. Soil Biol. Biochem. 120, 105–125. https://doi.org/10.1016/j.soilbio.2018.01.030. Burt, R., 2004. Soil survey laboratory methods manual. USDA-NRCS. Soil survey investigations report no. 42, version 4.0. Capra, G.F., Ganga, A., Grilli, E., Vacca, S., Buondonno, A., 2015. A review on anthropogenic soils from a worldwide perspective. J. Soils Sediments 15, 1602–1618. https://doi.org/10.1007/s11368-015-1110-x. Cerrejón, 2014. Informe de Sostenibilidad http://www.cerrejon.com/wp-content/2014min/IS%202014-min.pdf (accessed 02 january 2018). Chaves, L.F., 2010. An entomologist guide to demystify pseudoreplication: data analysis of field studies with design constraints. J. Med. Entomol. 47, 291–298. https://doi. org/10.1093/jmedent/47.1.291. Chessel, D., Dufour, A.B., Thioulouse, J., 2004. The ade4 package – I: one-table methods. R News 4, 5–10. https://doi.org/10.2307/3780087. Ciarkowska, K., Gargiulo, L., Mele, G., 2016. Natural restoration of soils on mine heaps with similar technogenic parent material: a case study of long-term soil evolution in Silesian-Krakow Upland Poland. Geoderma 261, 141–150. https://doi.org/10.1016/ j.geoderma.2015.07.018Get. Cooke, J.A., Johnson, M.S., 2002. Ecological restoration of land with particular reference to the mining of metals and industrial minerals: a review of theory and practice. Environ. Rev. 10, 41–71. https://doi.org/10.1139/A01-014. Cunha, L., Brown, G.G., Stanton, D.W.G., da Silva, E., Hansel, F.A., Jorge, G., Mckey, D., Vidal-Torrado, P., Macedo, R.S., Velasquez, E., James, S.W., Lavelle, P., Kille, P., de Indio, P., 2016. Soil animals and pedogenesis: the role of earthworms in anthropogenic soils. Soil Sci. 181, 1–16. https://doi.org/10.1097/SS.0000000000000144. De Paul Obade, V., Lal, R., 2016. Towards a standard technique for soil quality assessment. Geoderma 265, 96–102. https://doi.org/10.1016/j.geoderma.2015.11.023. Doblas-Miranda, E., Wardle, D.A., Peltzer, D.A., Yeates, G.W., 2008. Changes in the community structure and diversity of soil invertebrates across the Franz Josef Glacier chronosequence. Soil Biol. Biochem. 40, 1069–1081. https://doi.org/10.1016/j. soilbio.2007.11.026. Dolédec, S., Chessel, D., 1994. Co-Inertia analysis: an alternative method for studying species-environment relationships. Freshw. Biol. 31, 277–294. https://doi.org/10. 1111/j.1365-2427.1994.tb01741.x. Domínguez-Haydar, Y., Armbrecht, I., 2011. Response of ants and their seed removal in restoration areas and forests at El Cerrejón coal mine in Colombia. Restor. Ecol. 19, 178–184. https://doi.org/10.1111/j.1526-100X.2010.00735.x. Domínguez-Haydar, Y., Castañeda, C., Rodríguez-Ochoa, R., Jiménez, J.J., 2018. Assessment of soil fauna footprints at a rehabilitated coal mine using micromorphology and near infrared spectroscopy (NIRS). Geoderma 313, 135–145. https://doi.org/10.1016/j.geoderma.2017.10.032. Doran, J.W., Parkin, T.B., 1994. Defining and assessing soil quality. In: Doran, J.W., Coleman, D.C., Bezdicek, D.F. (Eds.), Defining Soil Quality for a Sustainable Environment. SSS Special Publication No 35. Soil Science Society of America, Madison, pp. 3–21. Dray, S., Chessel, D., Thioulouse, J., 2003. Co-inertia analysis and the linking of ecological data tables. Ecology 84, 3078–3308. https://doi.org/10.1890/03-01789. Escoufier, Y., 1973. Le traitement des variables vectorielles. Biometrics 29, 750–760. Frouz, J., Pižl, V., Tajovský, K., Starý, J., Holec, M., Materna, J., 2013a. Soil macro- and mesofauna succession in post-mining sites and other disturbed areas. In: Frouz, J. (Ed.), Soil Biota and Ecosystem Development in Post Mining Sites. CRC Press, Boca Ratón, USA, pp. 217–235. https://doi.org/10.1201/b15502-7. Frouz, J., Jílková, V., Cajthaml, T., Pižl, V., Tajovský, K., Háněl, L., Burešová, A., Šimáčková, H., Kolaříková, K., Franklin, J., Nawrot, J., Groninger, J.W., Stahl, P.D., 2013. Soil biota in post-mining sites along a climatic gradient in the USA: Simple communities in shortgrass prairie recover faster than complex communities in tallgrass prairie and forest. Soil Biol. Biochem. 67, 212–225. https://doi.org/10.1016/j. soilbio.2013.08.025. Granados, A.M., Barrera, J.I., 2007. Efecto de la aplicación de biosólidos sobre el repoblamiento de la macrofauna edáfica en la cantera Soratama. D.C. Univ. Sci., Bogotá, pp. 73–84. https://doi.org/10.1017/CBO9781107415324.004. [in Spanish]. Grimaldi, M., Oszwald, J., Dolédec, S., Hurtado, M.P., de Souza, Miranda I., Arnauld de Sartre, X., de Assis, W.S., Castañeda, E., Desjardins, T., Dubs, F., Guevara, E., Gond, V., Lima, T.T.S., Marichal, R., Michelotti, F., Mitja, D., Noronha, N.C., Delgado

192

Ecological Indicators 103 (2019) 182–193

Y. Domínguez-Haydar, et al.

Wanner, M., Dunger, W., 2002. Primary immigration and succession of soil organisms on reclaimed opencast coal mining areas in eastern Germany. Eur. J. Soil Biol. 38, 137–143. https://doi.org/10.1016/S1164-5563(02)01135-4. Welham, S.J., Cullis, B.R., Kenward, M.G., Thompson, R., 2007. A comparison of mixed model splines for curve fitting. Aust. New Zealand J, Stat. 49 (1), 1–23. https://doi. org/10.1111/j.1467-842X.2006.00454.x. Wilkinson, M.T., Richards, P.J., Humphreys, G.S., 2009. Breaking ground: pedological, geological, and ecological implications of soil bioturbation. Earth-Sci. Rev. 97, 257–272. https://doi.org/10.1016/j.earscirev.2009.09.005. Zangerlé, A., Hissler, C., McKey, D., Lavelle, P., 2016. Using near infrared spectroscopy (NIRS) to identify the contribution of earthworms to soil macroaggregation in field conditions. Appl. Soil Ecol. 104, 138–147. https://doi.org/10.1016/j.apsoil.2015.09. 014.

2012.01.008. Velásquez, E., Lavelle, P., Andrade, M., 2007a. GISQ, a multifunctional indicator of soil quality. Soil Biol. Biochem. 39, 3066–3080. https://doi.org/10.1016/j.soilbio.2007. 06.013. Velásquez, E., Pelosi, C., Brunet, D., Grimaldi, M., Martins, M., Rendeiro, A.C., Barrios, E., Lavelle, P., 2007b. This ped is my ped: visual separation and near infrared spectra allow determination of the origins of soil macroaggregates. Pedobiologia 51, 75. https://doi.org/10.1016/j.pedobi.2007.01.002 87. Walkley, A., Black, I.A., 1934. An examination of the Degtjareff method for determining soil organic matter, and proposed modification of the chromic acid titration method. Soil Sci. 37, 29–38. https://doi.org/10.1097/00010694-193401000-00003. Wang, Y., 2011. Smoothing Splines: Methods and Applications. CRC Press, Boca Ratón, pp. 354.

193