Assessing the Long-Term Ecosystem Productivity Benefits and Potential Impacts of Forests Re-Established on a Mine Tailings Site

Restoring sites disturbed by industrial activity to a forested condition can ensure the continued provision of economic and ecosystem services from these areas. Impounded mine tailings are particularly challenging sites, and positive benefits of establishing trees must be balanced against risks associated with metal contamination, ongoing tailings stability, and the possibility of acid mine drainage. We used a hybrid biometric modelling approach based on dendrochronological reconstruction to retrospectively (1980–2015) quantify productivity and carbon dynamics of pine plantations growing on impounded mine tailings at the Vale waste management facility near Sudbury, Canada. Historical reclamation practices had remediated conditions sufficiently to allow conifer plantation establishment in the late 1970s. The revegetated sites were highly productive, when compared to reference conditions based on site index, wood volume growth, and ecosystem production, congruent with other studies showing that forests on revegetated post mining sites can be highly productive. However, metal concentrations in the forest floor were high, and further research is warranted to evaluate ecosystem impacts. Due to the requirement for energy-intensive inputs, we estimated that it took 12 years or more to recover the emissions associated with the revegetation process through C accumulated in biomass and soil at the revegetated sites.


Introduction
Canada's forested regions contain significant mineral and energy resources, the extraction of which contributes to the fulfillment of society's needs for goods and services and makes a significant contribution to the national economy. Large amounts of minable resources have been and will continue to be extracted. One of the consequences of this is that the areas on which this activity takes place are subject to severe disturbance, which will require reclamation or restoration in order to be returned to a productive forested state and ensure the continuation of social license [1,2]. Returning industrially disturbed sites to a forested condition is critical to ensuring that these sites continue to provide a variety of both economic and ecosystem services. However, successful revegetation in boreal forests can take decades [3] and must take into account that forests are complex ecosystems consisting of long-lived plants where successional changes, through multiple pathways, can often take decades or centuries [4]. Revegetation must also be resilient to future changes in climate that will influence many ecosystem processes in Canadian boreal forests [5]. Activities undertaken in the forest sector are one of a potential suite of actions that could be undertaken to mitigate the potential effects of climate change [6]. The rehabilitation of industrially disturbed forest lands has been recognized as one possible activity to achieve these goals, but the potential of forests established on revegetated industrial sites to contribute has not been adequately quantified [7].
Impounded mine tailings are an example of an industrially disturbed site that presents revegetation challenges. Such sites are typically carbon (C) and nutrient depleted, have a massive soil structure, contain high concentrations of metal contaminants and can be acid-generating [8]. An example of such a site, which also has a long and rich history of revegetation trials, can be found on the Vale waste management facility within the city of Greater Sudbury, Ontario, Canada [9][10][11][12][13][14]. In the region surrounding the city of Greater Sudbury, smelting activities in the 20th century created large areas of metal-contaminated, phytotoxic land, which have in many cases been successfully revegetated [15,16], a process that continues to be implemented through a regreening program. The long-term effects of this program on ecosystem recovery continue to be evaluated and improved [17,18]. Similar to industrially damaged forests in the region, the historical reclamation practices on these impounded mine tailings consisted mainly of adding inorganic fertilizers, lime and/or mulch as needed to establish pioneer plant communities of grass and legume [10,11]. By the mid-1970s, conditions had been sufficiently remediated in some parts of the Vale waste management facility to enable a tree planting program to establish conifer plantations [10,11]. Additional considerations with any activity to revegetate tailings in humid and cold climates are maintenance of the stability of the tailings, prevention of acid mine drainage, and minimizing the transfer of metal contaminants to the surrounding environment [19,20].
In this paper, we evaluate wood volume production, accumulated ecosystem C stock in both biomass and soil, and annual rates of ecosystem production at 35-to 40-year-old jack pine (Pinus banksiana Lamb.) plantations established at the Vale waste management facility. Ecosystem C sequestration was examined in terms of three key ecosystem C fluxes: net primary production (NPP), net ecosystem production (NEP), and heterotrophic respiration (Rh). We used a hybrid biometric-modelling approach that relies on dendrochronological stand reconstruction data [21][22][23] to generate inputs to an ecosystem model, the Carbon Budget Model of the Canadian Forest Sector (CBM-CFS3, [24]) to perform these evaluations. We also compared the restored plantation sites to various validation and reference conditions. These included, (1) control sites at similar pine plantations located in managed forests in the Sudbury region, (2) average expected growth trajectories for jack pine forests in the province of Ontario as used in forestry management planning, (3) annual net ecosystem production estimates for jack pine forests regenerating after fire and harvest disturbance, as determined by eddy covariance studies, and (4) measurements of total soil organic C stock at the revegetated sites. The main objective of the study was to demonstrate that the productivity of restored forest stands, where typically few sites are available for generating estimates, can successfully be recovered using tree-ring approaches, and to determine how these compare to reference conditions. An additional objective was to provide insight into the distribution of metal contaminants (mainly Cu and Ni) in the revegetated tailings profile, and to determine if there was any evidence that forest restoration activity was generating risks to the stability of the tailings or increasing the possibility of acid mine drainage. Finally, we also comment on the question of whether there is any reasonable potential that forest restoration activities on highly disturbed industrial lands could result in any climate mitigation benefit when considering the GHG emissions of the amendments and the management activities needed to undertake the revegetation.

Study Region
Our data are from four study sites in the area of Sudbury, Ontario, Canada ( Figure 1). The revegetated sites are within the Vale waste management facility, and are 36-year-old jack pine plantations characterized by an organic layer that has accumulated over the years at the surface of oxidized, acid-generating Cu-Ni mine tailings. The tailings' surface is relatively flat, lying at an average elevation of about 320 MASL and covering an area of approximately 2500 ha. About half of the total area (1300 ha) is currently inactive (i.e., no longer receiving deposits of fresh tailings) and undergoing remediation. The total depth of tailings varies within the inactive impoundment from less than 1 m up to dozens of metres deep, depending on location. The water table varies seasonally, with lower-lying areas of the tailings being fully saturated during spring. In summer, the water table typically drops to about 2 m below the surface. At higher elevations (e.g., localized mounds of tailings), the upper, oxidized strata of the tailings remain unsaturated throughout the year.
Forests 2018, 9, x FOR PEER REVIEW  3 of 22 oxidized, acid-generating Cu-Ni mine tailings. The tailings' surface is relatively flat, lying at an average elevation of about 320 MASL and covering an area of approximately 2500 ha. About half of the total area (1300 ha) is currently inactive (i.e., no longer receiving deposits of fresh tailings) and undergoing remediation. The total depth of tailings varies within the inactive impoundment from less than 1 m up to dozens of metres deep, depending on location. The water table varies seasonally, with lower-lying areas of the tailings being fully saturated during spring. In summer, the water table typically drops to about 2 m below the surface. At higher elevations (e.g., localized mounds of tailings), the upper, oxidized strata of the tailings remain unsaturated throughout the year. The inactive tailings are subdivided for management purposes into areas designated by capital letters [9][10][11][12]. The study sites, which are located in an area referred to in previous publications [9][10][11][12] as the CD area, were established as part of a related project examining the biogeochemical evolution of acidic tailings under different reclamation strategies [20,25,26]. The CD area is about 132 ha in size including a central lake [9][10][11][12]. The two sampling sites were visually representative of the types of forested conditions currently found within the CD area, but they were not in a strict sense randomly located. Vegetation was re-established on these sites using the historical approach of applying fertilizers, limestone, and mulch in combination with a transition period with grass or legume crops prior to the establishment of a conifer plantation [9][10][11][12]. The details of the specific activities and amendments applied at the study sites, reconstructed to the best of our ability from historical reports [9][10][11][12], are summarized in Table S1. The main goal when revegetation activity initially began on this site was simply to establish vegetation to control blowing dust [11]. As a control, two additional sites were sampled in 44-and 45-year-old pine plantations, established following forest harvest in the managed forest that surrounds the Sudbury area, ~40-50 km away from the tailings site. A summary description of the basic characteristics of each revegetated and control site is provided in Table 1. The study region is located in the Boreal Shield ecozone of Canada's national terrestrial ecosystem classification system [27], and the sites are within ecoregions 4E and 5E of the ecological land classification system of the province of Ontario [28]. Over the 30 years prior to sampling, the mean annual temperature in the general study region was 4.5 °C (range 2.8 to 6.2 °C) with an average annual precipitation of 895 mm (range 685 to 1086 mm). The inactive tailings are subdivided for management purposes into areas designated by capital letters [9][10][11][12]. The study sites, which are located in an area referred to in previous publications [9][10][11][12] as the CD area, were established as part of a related project examining the biogeochemical evolution of acidic tailings under different reclamation strategies [20,25,26]. The CD area is about 132 ha in size including a central lake [9][10][11][12]. The two sampling sites were visually representative of the types of forested conditions currently found within the CD area, but they were not in a strict sense randomly located. Vegetation was re-established on these sites using the historical approach of applying fertilizers, limestone, and mulch in combination with a transition period with grass or legume crops prior to the establishment of a conifer plantation [9][10][11][12]. The details of the specific activities and amendments applied at the study sites, reconstructed to the best of our ability from historical reports [9][10][11][12], are summarized in Table S1. The main goal when revegetation activity initially began on this site was simply to establish vegetation to control blowing dust [11]. As a control, two additional sites were sampled in 44-and 45-year-old pine plantations, established following forest harvest in the managed forest that surrounds the Sudbury area,~40-50 km away from the tailings site. A summary description of the basic characteristics of each revegetated and control site is provided in Table 1. The study region is located in the Boreal Shield ecozone of Canada's national terrestrial ecosystem classification system [27], and the sites are within ecoregions 4E and 5E of the ecological land classification system of the province of Ontario [28]. Over the 30 years prior to sampling, the mean annual temperature in the general study region was 4.5 • C (range 2.8 to 6.2 • C) with an average annual precipitation of 895 mm (range 685 to 1086 mm).

Site Establishment
There is some uncertainty about when trees were established at the revegetated sites. One study mentions that conifers were first planted in the CD area of the Vale tailings in 1972, with a range of 5-15 thousand trees planted annually in the late 1970's and early 1980's [10]. Another suggests an even earlier date, in or around 1961 [12], likely coinciding with overall revegetation program inception rather than tree establishment at these locations. The mean breast height ring count in 2015 for both revegetated sites was 30 years, implying that trees reached 1.3 meters in height, on average, in 1985. A typical number of years for jack pine to reach breast height in this region is 8 years [30], giving an age in 2015 of 38 years and an establishment date of 1977. Ring counts obtained on samples very close to the ground contained a maximum of 33 rings, implying that trees were established in 1982. We used 1980 as the establishment date, which is between these two estimates, and is consistent with the timing of tree establishment at this site implied in [9]. Tree-ring counts have not previously been obtained at these sites, so we believe our establishment date estimate is close to correct.

Site Conditions
We provide here some additional description about the likely properties of the tailings before jack pine plantations were established. In the absence of original data, we relied on information reported in previous studies to determine the chemical and mineralogical composition of the tailings at the revegetated sites and to establish the presence of an oxidized tailings layer prior to tree establishment. The composition of tailings will vary with the properties of the mined ore. The average mineralogical composition of the Vale tailings was characterized by pyrrhotite (Fe 7 S 8 ) as the main sulfide mineral (~5.6 wt. % = 2.0 wt. %S), with pentlandite ((Fe, Ni) 9 S 8 , 0.5 wt. %) and chalcopyrite (CuFeS 2 , 0.3 wt.%) in lesser amounts [10], and various gangue silicates including quartz, albite, clinochlore, phlogopite and amphibole [20]. Given that chalcopyrite is the unique source of Cu, this would represent an average total concentration of 1039 mg Cu kg −1 in the tailings during that time. The CD area had a reported a sulfur content of 5.4 wt. % in the unweathered tailings (2 to 3 m below the surface) [31]. These data are in the same range as concentrations measured in unoxidized Vale tailings from an adjacent, more recent impoundment, with 1635 to 1785 mg Cu kg −1 , 2.35 to 5.35 wt. % of total S, and pH ranging from 3.8 to 5.0 [25]. Early grass seeding trials at this site failed mainly due to low pH, which required the addition of limestone [10,11]. This indicates that the pH of the tailings was initially very low, thus suggesting the presence of an oxidized layer. The occurrence of an oxidized layer at the surface is important to consider because the mineralogy of oxidized tailings (dominated by oxides and sulfates) differs significantly from that of the original tailings (dominated by sulfides). In Vale, oxidized tailings iron oxyhydroxide phases are abundant, with gypsum (CaSO 4 ·2H 2 O), goethite (α-FeOOH) and jarosite (KFe 3 (SO 4 ) 2 OH 6 ) as the main secondary minerals resulting from sulfide oxidation [20,32].

Tree Data
Data on past growth trajectories were obtained by dendrochronological stand reconstruction [21][22][23]. This involves obtaining increment core samples from all live and dead trees at fixed area forest plots, and can reliably generate estimates of past growth for several decades into the past, depending on dead tree persistence [33]. At the revegetated site, two 225 m 2 (15 m × 15 m) plots were sampled in October 2015. Identical plots were sampled at each control site in July 2016. A full census was made of all live and dead trees present at sampling. Height was measured on living trees; breast height diameter (DBH) on all trees. Two randomly oriented breast height cores were extracted from all living trees and a cross-sectional disc cut from all standing dead trees and downed logs. Samples were processed using standard dendrochronological methods, measured with WinDendro (Regent Instruments, Quebec, Canada) and cross-dated against site master chronologies developed from a subset of dominant trees at each site. All sample processing and measurement was conducted at the Northern Forestry Centre of Natural Resources Canada, Edmonton, Alberta. Due to decomposition, ring-widths could not be measured on 1.8% of trees at the control sites, for which ring-widths and year of death were estimated from the 3 larger and 3 smaller trees nearest in diameter. The mean absolute difference between true and estimated year of death using this method is 3.3 (SD 2.4) years when applied to jack pine trees with known year of death [34]. No trees were excessively decomposed at the tailings site. The timing of mortality and downed log persistence suggest 30-50 years into the past as a reliable reconstruction period [22,33]. Because our sites are plantations and generally younger than this persistence period, we considered our reconstruction of growth in this case to be reliable since stand establishment.

Soil Data
Soil samples were collected at both revegetated sites in order to measure the main metal contaminants (Cu, Ni, Fe) and the organic C content of the revegetate tailings profile. Four random locations, 20 m apart, were sampled at each site. At each location, a composite sample of the forest floor and surface litter was collected, then two 50-cm plexiglass columns were used to retrieve intact vertical core samples that included the barely decomposed litter and humus portions of the forest floor (FF) and the underlying oxidized tailings. These samples were transported on ice packs to the laboratory, where the cores were split into distinct layers: forest floor, the interface between forest floor and oxidized tailings, and the 0 to 5 and 5 to 15-cm layers of the oxidized tailings. No sample could be obtained from the bottom layer of unoxidized tailings due to tailing oxidation depth down to 80 cm depending on the locations and the occurrence of a hard pan layer at <20 cm depth. Total C in the samples was determined by combustion and infrared detection of CO 2 using an Eltra CS-2000 analyzer. The inorganic C was determined the same way after removal of organic C by combustion at 400-450 • C for 12 h. The organic C was obtained by difference. The pH was measured in distilled water with a 1:2 ratio (solids/water; w/v) for tailings and a 1:10 ratio for organic layers. Total Cu, Ni and Fe contents of tailings were determined for the CD-slope site only and were obtained by an acid digestion with HCl, HNO 3 , HF and HClO 4 [35] followed by ICP-AES.

Tree Diameter and Height
The annual DBH height (inside bark D ib ) was determined for each tree from the arithmetic mean of the cumulative ring-widths on each of two measured radii, correcting first where needed for the number of missing rings, distance to missed pith, and shrinkage. Diameter inside bark (Dib, cm) was converted to diameter outside-bark (D ob , cm) using: where k was assumed to be 0.964 for jack pine [22] and 0.951 for red pine [36]. An initial estimate of height (H i , m) was predicted from D ob (=D) using the Richards function where parameters a, b, and c were obtained from [34]. A correction factor derived from the ratio of the predicted and observed height at the time of sampling was then applied to the past estimates of height, such that the final estimated height was equal to measured height at the time of sampling. These corrected past heights (H c ) were then used to estimate past stand heights (SHT, m), and used to generate a second estimate of height using a model that also relies on stand height in addition to stand basal area (BA, m 2 ha −1 ), and stand density (TPH, stems ha −1 ) [37] where parameters θ, δ, β, ϕ and γ were obtained from [37]. This second estimate of height was also adjusted by the same tree-specific correction factor to generate a final estimate of height (H f , m), which reduces the error associated with estimating whole tree growth from increment cores taken at breast height only [38].

Wood Volume
Total volume (V T , m 3 ) was estimated using Honer's [39] standard volume equation Parameters C 1 , C 2 , and B 2 were obtained from [39]. Merchantable volume (V M , m 3 ) was determined using the adjusted squared diameter method [39] using the merchantability criteria assumed for Ontario by the CBM-CFS3 [24,40,41]. Tree values were summed to obtain stand-level volume (m 3 ha −1 ) and differenced to obtain volume increment (m 3 ha −1 year −1 ).

Carbon Budget Model of the Canadian Forest Sector (CBM-CFS3)
Reconstructed estimates of merchantable volume increment were used as inputs to the Carbon Budget Model of the Canadian Forest Sector (CBM-CFS3 [24]). The CBM-CFS3 uses these data in combination with other information to derive annual resolution estimates of past NPP, NEP, and Rh, as well as ecosystem C stock [22,23]. The model is described in existing publications [24,41]. In brief, it tracks three above-ground (foliage, merchantable stemwood, other biomass) and two below-ground (coarse roots and fine roots) biomass C pools by species group (softwoods and hardwoods). Models stratified by species and terrestrial ecozone are used to estimate above-ground biomass C from merchantable volume [40]. Below-ground biomass C is estimated from above-ground biomass C and species group [42]. Annual turnover of each biomass component is assigned to dead organic matter (DOM; standing and downed dead wood) C pools using litterfall transfer parameters (% year −1 ). Soil and DOM C stocks are tracked using eleven pools. Temperature-dependent annual decay decomposes the C in these pools, some of which is released to the atmosphere and the remainder transferred to a stable, slowly decomposing pool from which decay releases all C to the atmosphere.

Mean Annual Temperature
The CBM-CFS3 uses historical mean annual temperature to estimate past inter-annual variation in Rh. This was obtained from station records in the Sudbury area extracted from Meteorological Service of Canada's online archive using the Canadian Climate Data Scraping Tool [43].

Uncertainty Analysis
An uncertainty analysis was conducted in order to determine confidence intervals for the stand reconstruction estimates of ecosystem C stocks and fluxes [44]. Confidence intervals were derived from the 2.5th and 97.5th percentiles of n = 100 Monte Carlo simulations for each plot, where the simulations account for uncertainty in net biomass increment, litterfall, and soil and DOM C modelling parameters. Briefly, a triangularly distributed ±50% multiplier was applied to modelled net biomass increment to estimate uncertainty attributable to allometric models used to determine the input merchantable volume time series, and in the CBM-CFS3 to estimate biomass. Uncertainty attributable to biomass turnover and soil and dead organic C modelling parameters was estimated by varying 32 of these parameters, using probability distributions derived from literature values [44].

Soil and DOM C Initialization
An issue with modelling C dynamics on revegetated sites is the uncertainty about the initial soil and DOM C stock. These pools are usually initialized by repeated iterations of growth and stand-replacing disturbance until the difference between the total C in slowly decaying pools in successive iterations is less than 1%, followed by one additional disturbance reflecting the last known disturbance type. This procedure is assumed on average to yield soil and DOM C conditions consistent with the disturbance history of the stand [24] and is valid for the control sites in the managed forest, which were known to be established following harvest of a previous stand dominated by natural fire disturbance cycles. It is not valid for the revegetated sites, since they were not subject to this regime. Alternative procedures used in similar situations such as regeneration of forest on abandoned agricultural land derive these values from agricultural soil C databases [45] or use averages for nonforest soil types [41]. It would be reasonable to assume that the tailings were C-depleted at site establishment. However, because tree establishment was preceded by a period where various amendments were applied in combination with grass or legume crops, it cannot be assumed that soil and DOM C were entirely absent. To assess this key uncertainty, we conducted three sensitivity analysis simulations for the two revegetated sites, (1) zero initial soil C, (2) average initial soil C used by the CBM-CFS3 for afforestation sites [41], and (3) initializing soil C using the default procedure for natural stands [24].

Reference Estimates
We compared observed site index at each study plot against standard curves for jack pine used regionally in forest management [29]. We also compared observed wood volume yield to two standard yield curves used in forest management in Ontario, (1) "Plonski" [46,47], and (2) "Penner" [48]. We also compared our estimated NEP to independent estimates for regenerating jack pine stands at the Boreal Ecosystem Research and Monitoring (BERMS) study sites in Saskatchewan, Canada. In total, n = 6 eddy covariance sites (HJP75, HJP94, HJP02, F77, F89, and F98) originating in various years following harvest (H) or fire (F) were monitored from 2001 to 2005. Annual NEP values were obtained from tabular data appearing in publications detailing these results [49][50][51][52].
2.6.2. Field Measured Soil C Soil C stocks estimated from the samples taken at each revegetated site (Section 2.2) were compared against modelled output. Bulk density was not measured. For forest floor, we assumed 0.16 g cm −3 (SD 0.33), the average of the LFH layers for Brunisols in Canadian forest soils [53]. For the Vale oxidized tailings, where 73% of the mineral matter was comprised of fine clay, clay and silt, a lab estimate of 1.3 g cm −3 was obtained. However, this material would be more compacted in the field so we used this as the minimum, and 1.6 g cm −3 [54] as the maximum estimate for tailings bulk density. The average of forest floor and tailings was used as the estimate for the forest floor-tailings interface layer. The minimum and maximum bulk density estimates for tailings (95% CI in the case of the forest floor) were used to generate a range of estimates of total C for each layer and for the whole soil profile. The means and standard deviations from the n = 4 samples were used to calculate estimates of total soil C and associated confidence intervals for both the minimum and maximum bulk density estimates for each layer. The reported mean estimates are the average of the minimum and maximum estimates, and the reported confidence intervals are the widest possible range of values.

Tree Height and Site Index
The relationship between observed and predicted height was good when both DBH and stand characteristics (density, diameter, and stand height) were used as predictors (Figure 2). The CD site was an exception, and had trees that were on average 13.6% shorter than predicted. This site also had the lowest stand density and basal area. In general, the stem density of the revegetated sites was lower (<50%) of the density of the control sites. This was most visually evident by a much branchier growth form at CD than at the other sites. All sites were of broadly similar site quality. The estimated site index [29] was 18 and 23 m at a breast height age of 50 years for the control sites Cartier and Kukagami, respectively. It was 20 and 22 m for the revegetated CD site and CD-SLOPE sites.

Tree Height and Site Index
The relationship between observed and predicted height was good when both DBH and stand characteristics (density, diameter, and stand height) were used as predictors (Figure 2). The CD site was an exception, and had trees that were on average 13.6% shorter than predicted. This site also had the lowest stand density and basal area. In general, the stem density of the revegetated sites was lower (<50%) of the density of the control sites. This was most visually evident by a much branchier growth form at CD than at the other sites. All sites were of broadly similar site quality. The estimated site index [29] was 18 and 23 m at a breast height age of 50 years for the control sites Cartier and Kukagami, respectively. It was 20 and 22 m for the revegetated CD site and CD-SLOPE sites. Heights were predicted by models that use stand characteristics (density, basal area, and stand height) and DBH as predictors [37]. The straight line is a 1:1 line, and each study plot is identified by symbols and colours in the legend.

Wood Volume
Merchantable wood volume production, annual volume increment, and annual mortality for each study site are in Figure 3A-D. All sites are productive, with annual wood volume increment peaking at >10 m 3 ha −1 year −1 at the revegetated sites and Cartier and >15 m 3 ha −1 year −1 at Kukagami. Corresponding yield curves "Plonski" [46,47] and "Penner" [48] for fully stocked stands of site index 15, 20, and 25 m, which bracket the study plot estimates, are compared against observations in Figure  3E. Wood volume production was greater than "Penner" at all sites ( Figure 3E) and greater than "Plonski" at Kukagami. Wood volume production roughly tracks the site index 25 "Plonski" yield curve at the revegetated sites and Cartier. Heights were predicted by models that use stand characteristics (density, basal area, and stand height) and DBH as predictors [37]. The straight line is a 1:1 line, and each study plot is identified by symbols and colours in the legend.

Wood Volume
Merchantable wood volume production, annual volume increment, and annual mortality for each study site are in Figure 3A-D. All sites are productive, with annual wood volume increment peaking at >10 m 3 ha −1 year −1 at the revegetated sites and Cartier and >15 m 3 ha −1 year −1 at Kukagami. Corresponding yield curves "Plonski" [46,47] and "Penner" [48] for fully stocked stands of site index 15, 20, and 25 m, which bracket the study plot estimates, are compared against observations in Figure 3E. Wood volume production was greater than "Penner" at all sites ( Figure 3E) and greater than "Plonski" at Kukagami. Wood volume production roughly tracks the site index 25 "Plonski" yield curve at the revegetated sites and Cartier. ) and revegetated (CD (C), and CD-SLOPE (D)) plots sampled in this study. The legend in (C) applies to panels (A-D). Volume is in merchantable units, compiled to CBM-CFS3 standards, and plotted in panels (A-D) against calendar year on the x-axis. From ringcounts at ground height, the control plots were established in about 1971, the revegetated plots in about 1980. Total volume production (white squares) is plotted against the left y-axis. Annual volume increment (red circles) and mortality (grey bars) are plotted against the right y-axis. The volume increment time series in panels (A-D) are used as input in the CBM-CFS3 simulations to estimate ecosystem production. In panel (E), merchantable volume is plotted against stand age (years) on the x-axis for each study plot, relative to standard yield curves "Penner" ( [48], solid lines) and "Plonski" ( [46,47], dashed lines) for fully stocked jack pine stands with site index 15, 20, and 25 m, used regionally in forest management.

Ecosystem Production and Soil C Stock Estimates
Estimates of NPP, NEP, and Rh for the revegetated sites are in Figure 4A (CD) and Figure 4B (CD-SLOPE). Only results for the scenario assuming zero initial soil C are shown because estimated soil C stocks for this scenario best agreed with measured soil C ( Figure 5). Soil C estimates for the default and average scenarios were larger than measured ( Figure 5). Estimated NPP would be the same for all scenarios, but estimated NEP and Rh would vary. For scenarios with larger soil C estimates, Rh would be higher and therefore NEP lower because more C would be available for release by decay. In the last five years before sampling, NEP at CD averaged 250 g C m −2 year −1 1. At Figure 3. Wood volume and volume increment time series derived from stand reconstruction for the control (Cartier (A), and Kukagami (B)) and revegetated (CD (C), and CD-SLOPE (D)) plots sampled in this study. The legend in (C) applies to panels (A-D). Volume is in merchantable units, compiled to CBM-CFS3 standards, and plotted in panels (A-D) against calendar year on the x-axis. From ring-counts at ground height, the control plots were established in about 1971, the revegetated plots in about 1980. Total volume production (white squares) is plotted against the left y-axis. Annual volume increment (red circles) and mortality (grey bars) are plotted against the right y-axis. The volume increment time series in panels (A-D) are used as input in the CBM-CFS3 simulations to estimate ecosystem production. In panel (E), merchantable volume is plotted against stand age (years) on the x-axis for each study plot, relative to standard yield curves "Penner" ( [48], solid lines) and "Plonski" ( [46,47], dashed lines) for fully stocked jack pine stands with site index 15, 20, and 25 m, used regionally in forest management.

Ecosystem Production and Soil C Stock Estimates
Estimates of NPP, NEP, and Rh for the revegetated sites are in Figure 4A (CD) and Figure 4B (CD-SLOPE). Only results for the scenario assuming zero initial soil C are shown because estimated soil C stocks for this scenario best agreed with measured soil C ( Figure 5). Soil C estimates for the default and average scenarios were larger than measured ( Figure 5). Estimated NPP would be the same for all scenarios, but estimated NEP and Rh would vary. For scenarios with larger soil C estimates, Rh would be higher and therefore NEP lower because more C would be available for release by decay. In the last five years before sampling, NEP at CD averaged 250 g C m −2 year −1 1. At CD-SLOPE, it averaged 255 g C m −2 year −1 . The highest five-year average NEP at CD was 432 g C m −2 year −1 for the period 2000-2004 when the site was 21-25 years old. At CD-SLOPE, the highest five-year average was 387 g C m −2 year −1 for 1996-2000 when the site was 17-21 years old.

Ecosystem Production
Annual NEP is shown for the control sites since establishment in 1971 and 1972 in Figure 6A and for the two revegetated sites since 1980 in Figure 6B. Figure 6A,B also show annual NEP estimates from n = 6 eddy covariance sites regenerating after fire or harvest at BERMS jack pine sites in Saskatchewan [49][50][51][52], age matched to the same year that the revegetated or control sites would have been that age. Both the control site and eddy covariance estimates show NEP changing from a C source to a sink ~10 years post-disturbance. Before that, the control sites show larger C sources than eddy covariance. After, they show larger C sinks. Revegetated sites are a C sink almost immediately after establishment because Rh is very low (Figure 4). Thus, cumulative NEP is higher at the revegetated sites than the control sites at the same age ( Figure 6C). The revegetated sites are immediately a cumulative C sink, while control sites do not make up for the loss of C from the decay of biomass transferred to DOM following disturbance until ~20 years after establishment ( Figure 6C).

Ecosystem Production
Annual NEP is shown for the control sites since establishment in 1971 and 1972 in Figure 6A and for the two revegetated sites since 1980 in Figure 6B. Figure 6A,B also show annual NEP estimates from n = 6 eddy covariance sites regenerating after fire or harvest at BERMS jack pine sites in Saskatchewan [49][50][51][52], age matched to the same year that the revegetated or control sites would have been that age. Both the control site and eddy covariance estimates show NEP changing from a C source to a sink~10 years post-disturbance. Before that, the control sites show larger C sources than eddy covariance. After, they show larger C sinks. Revegetated sites are a C sink almost immediately after establishment because Rh is very low (Figure 4). Thus, cumulative NEP is higher at the revegetated sites than the control sites at the same age ( Figure 6C). The revegetated sites are immediately a cumulative C sink, while control sites do not make up for the loss of C from the decay of biomass transferred to DOM following disturbance until~20 years after establishment ( Figure 6C).  [49][50][51][52]. The stand age at the time of EC observations is matched to the same age of the control sites or of the revegetated sites (in this case, age equals years since reclamation).  [49][50][51][52]. The stand age at the time of EC observations is matched to the same age of the control sites or of the revegetated sites (in this case, age equals years since reclamation).

Ecosystem C Stock
Despite accumulating at a faster rate, the revegetated sites in 2015 had only~38% of the total ecosystem C of the control sites at that age ( Figure 7A). Biomass C stocks at the revegetated sites werẽ 80% of the control sites ( Figure 7B) and soil and DOM C stocks were~12% ( Figure 7C).

Ecosystem C Stock
Despite accumulating at a faster rate, the revegetated sites in 2015 had only ~38% of the total ecosystem C of the control sites at that age ( Figure 7A). Biomass C stocks at the revegetated sites were ~80% of the control sites ( Figure 7B) and soil and DOM C stocks were ~12% ( Figure 7C).   Tailings Profile   Table 2 reports pH, total organic C and concentrations of Cu, Ni, and Fe in the 0-15 cm layers above the cemented layer of the revegetated tailings at the CD-slope site. Relative to Sudbury regional soils [55], all samples were elevated in Cu, Ni and Fe, but pH of the revegetated tailings was similar to Sudbury soils (typically < 5.0). Cu and Ni are reported as the main contaminants of concern, and they are dominantly associated with goethite in the oxidized tailings [20,25]. Fe represents the Fe-oxyhydroxide matrix of the oxidized tailings. Unreclaimed, oxidized Vale tailings typically have very acidic pH and are depleted in organic C (pH 2.47 and TOC 1.77 g kg −1 ; Table 2). In the revegetated tailings, a marked improvement in pH and organic matter content was noted in the 0-15 cm layer of the oxidized tailings (pH > 5.5 and TOC > 4.3 g kg −1 , Table 2). The average Fe concentration was less than half of what would typically be found in unreclaimed, oxidized Vale tailings (10.9 vs. 24.1 wt. % Fe respectively; Table 2). In these oxidized tailings, goethite and jarosite are usually the dominant Fe-bearing phases but characterization by bulk X-ray diffraction did not detect either of these minerals, suggesting mineralogical evolution of the reclaimed system [20,25]. The average Cu and Ni concentrations in the 0-15 cm layer of the oxidized revegetated tailings were 324 mg Cu kg −1 and 313 mg Ni kg −1 , and these concentrations were in the same range as those measured at a control, unreclaimed tailings site from an adjacent, more recent impoundment (305 mg Cu kg −1 and 477 mg Ni kg −1 ; Table 2). However, Cu and Ni accumulated in the organic layer at the revegetated tailings sites to 2 times greater concentrations than in the oxidized tailings (776 mg Cu kg −1 , 805 mg Ni kg −1 ; Table 2). Table 2. Chemical properties of soil from the CD-Slope site including pH, total organic C (TOC), Cu, Ni and Fe concentrations (mean ± SE, n = 4) in the forest floor (FF) and oxidized tailings (ox). As a comparison, data from the surface layer of control, unreclaimed oxidized tailings (control ox 0-5 cm) and from Sudbury regional soils [55] are also reported. Only TOC was measured at CD site.

Discussion
We demonstrate how a hybrid biometric modelling approach based on dendrochronological stand reconstruction can be used to retrospectively quantify the productivity and C dynamics of pine plantations at revegetated sites growing on impounded mine tailings at the Vale waste management facility near Sudbury, Canada. Even though the revegetated sites were not established with the idea that they would develop into economically harvestable stands, they were highly productive when compared to reference conditions on the basis of site index, wood volume growth, and ecosystem production. These results are congruent with many other studies that show that forests on revegetated post mining sites can be highly productive [56,57]. For example, white pine plantations growing on reclaimed surface mine land in Virginia was projected to have a site index higher than natural white pine stands in that region [58], while three different boreal species in Alberta growing on restored oil sands mining sites have similar site indices to undisturbed natural stands [59]. Like previous studies [7], total C stock of the revegetated site had not yet reached the same levels as what was present in analogous managed forest stands growing under natural conditions (Figure 7). However, the cumulative NEP since stand establishment was larger for the revegetated site than for the analogous managed forest stands in the region ( Figure 6C). Annual NEP was also larger than estimates derived from eddy covariance studies in young disturbed forest stands elsewhere in Canada [49][50][51][52] (Figure 6B). The primary reason for this was low (or no) initial soil and DOM C stocks when the trees were established. Therefore, little or no soil or DOM C was available for decay, and soil C built up over time [7,60]. Unlike sites that were previously peatlands with very large soil C stocks [61], it may be possible in pine ecosystems to eventually reach similar levels of ecosystem C stock as present prior to industrial disturbance [62].
The accumulation of C in soils developing in pine ecosystems at post-mining sites in Poland has been estimated to range from 0.73 to 5.26 Mg ha −1 year −1 depending on the substrate [63]. This is larger than the average of 0.5 Mg ha −1 year −1 that has on average accumulated at the revegetated sites in this study, but our sites were also accumulating less biomass than in the Polish study [63], and would have been considered potentially less fertile due to the requirement to ameliorate pH by adding lime.
The modelling approach predicts that little C would be lost from the system through decay, and therefore both the annual C sequestration and total ecosystem C accumulation at the revegetated sites is higher than at control sites located at pine plantations established following forest harvest in the managed forest that surrounds the Sudbury area. Modelled estimates of soil C that assumed that the sites started in a very C depleted state were most consistent with measurements of soil C obtained at the time of sampling. The observed organic C content at the revegetated sites relative to observations made 20 years previously [12], particularly in the forest floor, are consistent with an increase in soil C stock over time. Similar to an analogous study of a post-mining site in the Czech Republic [64], most soil C accumulation was in situ, originating as input from biomass turnover and litterfall. The decomposition of this litterfall is an important ecological process contributing to soil restoration [65]. In general, it is apparent that soils make a significant contribution to total C accumulation at revegetated sites on land that was formerly degraded [60]. The close match between modelled and measured values suggests that the soil and DOM C dynamics at the restored sites, including decomposition and transformation rates, were in close approximation to what is assumed by the soil and DOM C dynamics sub-module of the CBM-CFS3, which is parameterized from natural forest conditions [24]. It also provides further evidence that the soil and DOM C stock and flux estimates in the model are highly dependent on biomass input [66]. A small amount of soil C was likely initially present, originating from the mulch, grass, or legume crops present during the transition period prior to the establishment of a pine plantation [10,11]. This could be inferred from the difference between the measured soil C at the time of sampling and the model predicted soil C assuming the initial amount was zero.

Evaluation of Risks to Tailings Stability
This study shows that acid-generating tailings reclaimed with jack pine were highly productive 35 to 40 years after planting. However, successful reclamation should also aim to minimize acid mine drainage (AMD) and metal contaminant transfer into the environment (e.g., vegetation uptake). While inputs of C and other nutrients are key to successfully establishing vegetation on tailings, their addition or accumulation in these mineral wastes may also change the cycling and stability of metal contaminants in ways that could be deleterious to the health of the surrounding ecosystem. For example, the concentrations of Cu (776 mg kg −1 ) and Ni (805 mg kg −1 ) in the forest floor at the reforested CD-slope site both exceed the Canadian soil quality guidelines for industrial sites (91 and 89 mg kg −1 for Cu and Ni, respectively [67,68]). As the guidelines are based on toxicity data from plants and invertebrates exposed to affected soils, the levels of Cu and Ni measured in the reforested CD tailings suggest the potential for metal transfer to the food chain that could, for example, have toxic effects on invertebrate species [69,70]. A previous study on metal uptake by vegetation at this site [12] reported that Cu and Ni concentrations were higher in grasses, tree foliar tissues and mosses growing on tailings than at control sites. For grasses, Cu and Ni concentrations were also higher in roots than shoots. However, differences in metal accumulation depended on plant species, and their study showed comparable Cu and Ni concentrations in jack pine needles established on tailings at the CD area compared to nontailings sites for samples taken in 1996 [12]. In light of these studies, additional research is warranted to evaluate the impact of elevated metal concentrations on the ecosystem at these sites. On-going studies are also clarifying the degree of microbial recovery using molecular biology techniques to assess microbial community composition and function.
In addition to potential ecosystem effects, the reforested tailings profile exhibited extensive oxidation (down to >80 cm). This indicates that the revegetation approach did not prevent the oxidation of the unweathered tailings and AMD, though it was not expected that it would have been fully successful at doing so when the revegetation activities were initiated [10,11]. In line with the absence of goethite and jarosite as dominant Fe-bearing secondary mineral phases in the reforested tailings, chemical analysis indicated a significant depletion of Fe in the oxidized tailings compared to the non-reforested tailings control site ( Table 2). These observations likely reflect on-going podzolization which is favored under cold humid climates and coniferous vegetation. The beginning of a podzolic layer was already observed at these sites prior to the establishment of the jack pine plantations [10,11]. As the tailings are colonized by vegetation, organic acids are produced by decomposition of litter and exudation of roots, fungi and microbes. These organic acids contribute to mineral weathering by complexing Al and Fe and promoting their migration downwards where they can precipitate and form a hard pan [71], which was observed at the revegetated sites below 20 cm depth in the oxidized tailings. A common driving force for podzolization under coniferous vegetation is the occurrence of ectomycorrhizal fungi, which are believed to be responsible for the upward translocation of nutrients and other elements from the weathering mineral layer to the plant roots [72]. This phenomenon, coupled with potential dust contamination from surrounding barren tailings and vegetation uptake, may have contributed over the years to the accumulation of Cu and Ni in the tailings forest floor. These results certainly suggest that there are environmental risks to this particular tailings revegetation strategy, even though the site is presently a productive forest. An alternative approach that is currently being tested at this site is the use of oxygen consuming organic covers [13,14,20]. This approach, while potentially accelerating the rehabilitation of the tailings, also requires more assessment of its suitability to prevent acid mine drainage particularly if more economical thin covers are used [20]. The efficiency of a reclamation strategy in minimizing AMD should especially be considered at abandoned acid-generating sites where mine drainage is not collected for treatment like it is at this study site.

Assessment of GHG Mitigation Benefits
We examined the question of whether there is any reasonable potential that the revegetation of these sites could provide any GHG mitigation benefit. The rehabilitation of industrially disturbed forest lands has been recognized as an activity that could provide such benefits, but these have not often been adequately quantified [7]. Such an assessment must also consider the GHG emissions of the amendments and management activities needed to undertake the revegetation, ideally through a formal life cycle analysis. A similar study of an oak plantation established at a former lignite mine in Italy showed that tree growth was sufficient to offset the GHG costs of plantation management four years after establishment [73]. However, this was not an entirely analogous situation as the site in that study was actively managed with periodic thinning and required only standard mechanized forestry operations to establish trees, and not intensive input of high GHG emission amendments as was required at this study site.
The previously mentioned uncertainty about the precise establishment date of these sites, which was inconsistently reported in previous studies relative to the actual ages of the trees measured in this study using dendrochronological methods, makes it difficult to determine the precise timing and nature of various activities undertaken at these sites. However, some general assumptions may be made based on available published information [9][10][11]. As previously described, the amendments and activities required before trees could successfully be established, included the spreading of agricultural limestone, straw mulch, seeding with grasses and legumes, and fertilizer application. Some information about typical rates of application can be determined. Average rates of 25 Mg ha −1 for agricultural limestone, and 750 kg ha −1 of fertilizer are suggested [11], but these are also not consistent in different publications [9,10]. Supplementary Table S1 provides an outline of the best available estimates of the activities required to revegetate these sites based on this published information [9][10][11], as well as the potential GHG emissions associated with these activities, obtained from the literature [74][75][76]. Accounting only for the GHG emissions of the activities described in Supplementary Table S1 gives an estimate that could range from 2.4 to 5.7 Mg CO 2 ha −1 for activities and amendments required to revegetate these sites. Based on the maximum value in the range (5.7 Mg CO 2 ha −1 ), the accumulated NEP at the revegetated sites recovers this C in 7 to 12 years, which is comparable to analogous estimates for forest bioenergy projects that use harvesting residues under some circumstances [77,78]. However, significant uncertainties in these estimates remain. They do not account for all of the activities required in the process, particularly activities indicated as occurring but for which no estimate of the amount or details of the material applied are given, for example, site leveling prior to revegetation or the application of straw mulch and chemical binders [10,11]. There are also key uncertainties in the amounts of amendments applied. For example, some publications give an estimate of the amount of fertilizer applied as 750 kg ha −1 in total [10,11], while another gives an estimate of an annual rate of 300 kg ha −1 , which could have begun as early as 1958 [9]. In addition, there is the possibility that the reason for the discrepancy in precise timing of site establishment is that several attempts were made to establish these sites, and previous attempts failed. Therefore, the estimated emissions associated with the revegetation activities should be considered minimums. The estimate also does not account for non-CO 2 emissions associated with fertilization, and depends on a baseline scenario that assumes the sites are not revegetated.

Conclusions
In this study, we examined the recovery of revegetated sites located on impounded mine tailings at the Vale waste management facility near Sudbury, Ontario, Canada. In terms of wood volume growth and annual NEP, the productivity at revegetated sites was similar to or greater than that from regional control sites and other reference conditions to which the sites could be compared, despite a planting density that is lower than what is typically used in forest management regionally. Less C had accumulated in the soil and DOM at the revegetated sites than at the regional control sites, and the amount that had accumulated was consistent with low initial soil C stocks. A hybrid-biometric modelling approach that incorporates tree-ring data into an ecosystem model (CBM-CFS3) was successfully used to recover these ecosystem C dynamics, and could be applied to other historically reclaimed industrial sites. Can the results with respect to forest productivity at the revegetated sites in this study be indicative of successful reclamation? Relative to the initially stated goal when ecosystem restoration activities began on this site, which was simply to control blowing dust, it could be said that the objectives were achieved. Further investigation is needed to determine if metal accumulation in the forest floor originates from metal uptake by the trees and whether these elevated concentrations are having environmental impacts. This would be particularly important if there was a desire to use the wood biomass in forest products, though such use in not currently contemplated. Finally, despite the large quantity of GHG-intensive inputs required, the revegetation of industrially disturbed sites in this study may possibly have climate change mitigation benefits with similar C debt repayment times as some forest residue-derived bioenergy projects. However, significant uncertainties in these calculations remain, and the amount of time required to recover the GHG emissions associated with the reclamation activities and amendments required to revegetate the sites in this study are more likely to be longer than what was calculated than they are to be shorter.