Paludification of boreal soils reduces wood decomposition rates and increases wood-based carbon storage

. Over long time periods, paludification reduces aboveground productivity resulting in forest retrogression. Paludified forests are typified by intense accumulation of the soil organic layer and a reduction in soil temperatures and nutrient availability. En route to paludification, early successional forests experience large inputs of deadwood biomass during the senescence of the post-fire cohort, much of which may be entombed in this rapidly growing soil organic layer. Here we examined the effects of paludification across a . 2000-year chronosequence of black spruce forests on wood decomposition using three complementary approaches. We (1) repeatedly measured wood density of logs through time, (2) utilize a time-series of logs that varied in time since death, and (3) estimate woody biomass at the stand level as it progresses from live trees to snags, logs and ultimately to buried or decomposed deadwood. Together these approaches demonstrated a 6– 7-year delay before the onset of rapid decomposition. We also found strong evidence that paludification results in a large proportion of logs becoming buried in the soil organic layer. Stand level modeling indicates that the rates of accumulation of buried deadwood were greatest following the senescence of the post-fire cohort when both soil organic layer build-up and creation of deadwood peaked. Following this period of high deadwood creation, stands enter a retrogressive state whereby productivity continues to decline albeit more slowly. Continued losses in woody carbon biomass from trees during this retrogressive state are offset by lower wood decomposition rates and a high biomass of accumulated buried deadwood, essentially stabilizing the wood based carbon budget in these ecosystems. We recommend that partial cutting be conducted prior to or near the senescence of the post-fire cohort to improve emulation of natural forest succession in terms of both live tree and deadwood biomass. Furthermore, deadwood during this period has an extremely short residence time and the dynamics of deadwood should recover much quicker than if harvesting is conducted later in succession when there is less live tree biomass and deadwood has longer residence times.


INTRODUCTION
In the absence of disturbance, forest ecosystems can enter a retrogressive stage where aboveground productivity declines (Wardle et al. 1997, 2004, Richardson et al. 2004, Vitousek 2006).In northern forests, declines in productivity can be caused by paludification, whereby a thick soil organic layer, primarily in the form of Sphagnum mosses, accumulates resulting in colder, wetter soils with reduced nutrient availability, which ultimately reduces tree growth rates (Simard et al. 2007).With time, paludified forests develop into forest peatlands.Peatlands in general constitute 3% of global land area and contains 600 gigatons of carbon, more than a third of the world's pool of soil carbon (Woodwell et al. 1989, Yu et al. 2011).However, most of these peatlands are the result of edaphic paludification, where site topography and poor drainage result in peatland formation.Successional paludification, which has received relatively less attention, is the result of forest succession and fire cycles.Forests prone to successional paludification require high intensity fires to reset succession; whereas other disturbances (i.e., low intensity fire, insects and wind) tend to accelerate the paludification process.It is difficult to estimate the area of forest prone to successional paludification across the boreal; however, as the fire return interval in many parts of the boreal continues to lengthen (Bergeron et al. 2004), forest retrogression as a result of successional paludification will become more prevalent.Paludification of forested ecosystems and associated accumulation of soil organic material will increase their importance for carbon storage reservoirs and has the potential to drastically change deadwood dynamics.
Deadwood dynamics, like other aspects of boreal forest succession, are dependent on the frequency and intensity of disturbance events, site productivity and woody decomposition rates (Harmon et al. 1986, Bergeron et al. 2004).Initial volumes of deadwood are determined by disturbance events that reinitiate succession.For example, Siitonen (2001) reported increasing volumes of deadwood following wildfire (from 110 m 3 /ha to over 400 m 3 /ha) in spruce-dominated stands.These volumes depended on the amount of deadwood consumed by the fire, and the mortality rate of living trees.As initial deadwood decays, deadwood volumes are relatively low, until self-thinning increases the rate of deadwood deposition.In boreal forests, selfthinning usually begins 75-100 years after the disturbance, depending on site characteristics (Harvey et al. 2002, Harper et al. 2005).With time, input and decay rates can stabilize resulting in an equilibrium volume, which is common in many old growth stands (Tyrrell andCrow 1994, Duvall andGrigal 1999).
While our understanding of deadwood dynamics and succession are progressing (Krankina and Harmon 1995, Siitonen 2001, Brais et al. 2005), recent studies have drawn attention to the importance of buried wood (Manies et al. 2005, Hagemann et al. 2009).Rates of deadwood burial are determined by amount of ground contact, organic layer depth, canopy cover, and diameter of the log (Stenbacka et al. 2010).High rates of burial in some forest types can result in large volumes of deadwood under the surface.For example, Manies et al. (2005) found that 8-20% of woody biomass was buried in black spruce forests in Manitoba.Hagemann et al. (2009) reports that in high-boreal black spruce forests of Labrador, Canada, buried deadwood volumes exceeded aboveground volumes by 50% in old growth stands and .400% in 34-35 year old post-harvested stands Furthermore, the cool wet conditions under the moss layer reduce decay rates to the point where buried deadwood may persist until consumed by high intensity fire (Hagemann et al. 2009).When included, this large and persistent mass of buried deadwood increases total carbon in forest carbon budgets (Moroni et al. 2010).
Here we studied the effects of paludification on the input and losses of deadwood by measuring changes in biomass of live trees and deadwood across a chronosequence spanning 2375 years of stand development.We used three complimentary approaches to estimate loss from wood decay.We (1) repeatedly measured wood density of individual logs through time, (2) utilized a time-series of logs that varied in time since death, and (3) estimated woody biomass at the stand level as it progresses from live trees to snags, logs and ultimately to buried or decomposed deadwood.We then used these stand-level models to predict long-term woody carbon storage in paludified black spruce stands.We hypothesized that decreased productivity from paludification (Simard et al. 2007) will equate to lower rates of deadwood input, and the cold, wet soils will lower rates of wood decay compared to younger less paludified forests.We further expected a significant proportion of deadwood will also be lost to burial in the organic layer and as this buried wood accumulates with time, it will represent a large reservoir of stored carbon.

Study sites
The clay-belt region of Ontario and Quebec, Canada, covers ;12.5 3 10 6 ha and forms the second largest peatland in the world (Gorham 1991).It is part of the precambrian shield and is composed primarily of clay deposits left by proglacial lakes Barlow and Ojibway (Vincent andHardy 1977, Veillette 1994).The cold climate and flat topography, combined with clay soils, makes this region prone to paludification (Lavoie et al. 2005).Mean annual temperatures in the region are 08C with average precipitation values of 897 mm (33% falling as snow; Environment Canada 2011).
Our study was conducted along a chronosequence of stands in the northern part of the claybelt (49800 0 -50800 0 N; 78830 0 -79830 0 W) originating from high severity fires and ranging from 60 years to 2360 years in age (Fig. 1; Lecomte et al. 2006b).This chronosequence has been wellstudied in terms aboveground live biomass (Lecomte et al. 2006b), understory vegetation (Lecomte et al. 2005), diversity of mosses (Fenton and Bergeron 2006), lichens (Boudreault et al. 2009), litter-dwelling arthropods (Paradis and Work 2011) and forest productivity (Simard et al. 2007).We added three additional stands (vis.Pui, Fen, Gau) to this chronosequence, which were the uncut controls from an experimental study of partial cutting in the area (Jacobs andWork 2012, Fenton et al. 2013).A total of 15 stands were selected to represent the greatest possible time between stem exclusion (34-96 years), senescence of the post-fire cohort and understory reinitiation (96-164 years) and old growth stages (.164 years;Harper et al. 2005).

Stand-level deadwood dynamics model
Our approach for determining inputs and losses of deadwood across the chronosequence was based on estimations of biomass of live trees, snags and fallen logs which either decay on the surface or are buried beneath Sphagnum mosses where decomposition rates are thought to be greatly reduced (Fig. 2).Biomass of live trees and deadwood was estimated as the product of wood volume and wood density.Transitions of woody biomass among different deadwood pools were then modeled as a function of time since fire (TSF) using estimates of tree growth (for live trees), tree death rates (for snags), snag fall rates and decomposition rates measured both aboveand belowground (for logs).

Estimation of live tree biomass across the chronosequence
To estimate live tree biomass across the chronosequence, we converted individual tree volumes to mass using empirically derived estimates of volume and woody density from the chronosequence.Live tree biomass was assessed using circular plots with a diameter of 11.28 m (400 m 2 ).We measured each stem over 3 cm diameter at breast height (DBH), and recorded the species and whether the stem was alive or dead.Height of live trees was determined as a function of tree diameter using a calibrated Chapman Richards taper function (Zhao-gang and Feng-ri 2003).Individual parameters for the taper function were based on 107 black spruce trees from nine different stands (Boudreault et al. 2009) and fit using model comparisons with Akaike's information criterion corrected for finite sample sizes (AIC c ). Volume of each individual live tree, in the circular plot, was then calculated as a truncated cone extending from the base of the tree to a height where stem diameter was 5 cm.This necessarily excluded branches and smaller parts of the tree from our estimates of live tree volume that decay at different rates than larger diameter pieces.To convert volume to biomass we modeled wood density.We used measured wood density of living trees from nine stands and fit these densities with generalized linear models (GLM) and AIC c to determine the relationship between tree diameter, time since fire (TSF) and wood density.Biomass was then calculated by multi-plying predictions from the best-fit GLM model for each tree sampled across the chronosequence by live tree volume.

Estimation of snag biomass across the chronosequence
To characterize stand-level density of standing deadwood (snags) within each plot, we established a rectangular sub-plot (20 3 10 m) and recorded the diameter of each snag and whether the top was broken.In cases where the breakage was below 5 cm diameter, we estimated the height in the field.Heights of intact snags and snags broken above 5 cm diameter were estimated using the same taper function used for live trees.Snag volume was likewise estimated as a truncated cone.Wood density of snags were considered to be equal to the live wood density, as other studies have not found a significant relationship between wood density and time since death in black spruce snags (Boulanger and Sirois 2006, Boulanger et al. 2011, Angers et al. 2012).

Estimation of log biomass across the chronosequence
To estimate total log biomass, we measured volume, decay stage, density and mass of downed deadwood (logs) in two star plots within each stand (Sta ˚hl et al. 2001).Each star plot consisted of three 20-m transects radiating from a common midpoint and separated by 1208.For each log .5 cm diameter intersecting a transect, we recorded (1) the decay class using a five class system (modified from Maser et al. 1979), (2) the percentage of moss covering the log as an average of two meters on each side of the intersection point, and (3) the diameter with a diameter tape or caliper.Logs were classified as buried if moss cover was greater than 50%.We then cut a 2-5 cm cross-section from the log using a chainsaw.We measured the minimum and maximum thickness of each cross-section with a calliper and calculated the volume based on a cylinder.We then dried the cross-sections at 658C to constant mass and calculated wood density by taking the dry mass divided by the volume.We calculated stand level volumes of downed deadwood using the Van Wagner (1968) formula for overall volume and for each 0.1 g/cm 3 density class.
In two of the stands (N50 and N6, 373 and 710 years since fire, respectively), we excavated a 10 m by 30 cm trench to the depth of the mineral soil and measured all pieces of deadwood that were missed in the aboveground transects.For each piece, we measured the diameter and depth from the soil surface.We attempted to sample the oldest stands (N16 and N20) using these methods but, due to accessibility issues, it was not possible.

Measurement of organic layer depth across the chronosequence
We measured the depth of the organic layer along one transect of the star plot conducted for the estimation of log biomass.At one meter intervals along the transect, we pushed a steel rod into the ground until significant resistance was met and recorded the depth.

Estimation of decay rates of logs
Estimation of decay rates was based on changes of wood density of logs with increasing time since death of the tree.We established a sequence of logs varying in time since death in each of eight stands.Stand ages, based on time of sampling (2008), were 60, 94, 133, 134, 183, 373, 718 and 2365 years since the last fire.When possible we sampled the wood density of three age classes of logs: Age classes of logs were (1) zero years since death (freshly cut trees we felled in 2008), (2) four years since death (trees felled during study of Boudreault et al. [2009]), and (3) 6-10 years since death (trees felled during the study of Lecomte et al. [2006b]).
We initially sampled four (94, 134, 373 and 718 years since fire) of the nine stands in July of 2008.We cut multiple cross-sections from each log and prepared 50 cm segments of the youngest and oldest log age classes for later remeasurement.The number of cross-sections for each log age class were as follows: (1) for logs in the youngest age class (0 years), we cut four cross-sections 50 cm apart, starting 50 cm from the base of the tree; (2) in the second age class (4 years), we cut two cross-sections 50 cm apart, starting 50 cm from the base of the tree; and (3) in the oldest age class (6-8 years), we cut between six and 10 crosssections.Some of the trees in this oldest age class had cross-sections previously removed at 1-m intervals at time of felling.For these trees we cut two cross-sections for each interval for five corresponding sections.When we could find trees still intact, we cut six cross-sections at 50cm intervals.
Preliminary analyses indicated that by using additional stands we would decrease the variability in decay estimates; therefore we sampled three additional stands (60, 133, 373 years) in September of the same year.As in the first collection trip we sampled the same three age classes; however, we cut just two cross-sections, 50 cm from the base of each tree separated by 50 cm.Further analysis revealed possible effects of stand age on decomposition rates, which motivated us to add two extremely old stands to the study (1595 and 2365 year since fire).These stands were not sampled by Boudreault et al. (2009), leaving us with only the oldest and youngest age classes for logs.Furthermore, we were unable to find the trees cut by Lecomte et al. (2006b) in the 1595-year-old stand, leaving us with only a measurement of the youngest age class.Although the addition of these sites prohibit a fully balanced statistical comparisons, these data for living wood density, as well as the logs in the oldest stand, add to the overall value of this study.We compensated for these differences within log age classes and between stands by using mixed-effects models which are more robust to unbalanced sampling designs (Pinheiro and Bates 2000).
After measuring wood density for all crosssections, we sealed both ends of 50-cm segments from the youngest and oldest log age class with paraffin wax to reduce the rate of drying, which simulates the condition of longer logs.Logs were placed on the forest floor in the stand where they were sampled.Two years following initial sampling, we remeasured these log segments, for wood density, moss cover and canopy openness over each one using a Model A spherical densiometer (Lemmon 1956).

Data analysis
We compared changes in wood density of the same logs between years with repeated measures ANOVA (RM-ANOVA), using the aov function in R 2.14.1 (R Development Core Team 2011).Aberrant logs, where measured wood density increased, were removed from the analysis.We analyzed moss growth on these logs using LME with log nested in stand as the random effect and log diameter, stand age and canopy openness as the fixed effects.Interactions were not examined and models were compared with AIC c .
We initially modeled changes in our estimation of biomass of live trees, snags and logs, as well as changes in the depth of the organic layer using linear regression with the log (TSF) as the explicative variable.We then tested these models for the presence of a non-constant regression parameter (i.e., break-point) using the Davies test (Davies 1987).When a break-point in the regression was detected, we used segmented linear regression using the function segmented in the segmented library (Muggeo 2008) in R.This analysis estimates unknown break-points (Muggeo 2003) and corresponding regression coefficients.We removed the youngest stand from this analysis, as it was the only stand that was likely still gaining biomass at the time of sampling.
We estimated decay rates of logs by modeling changes in wood density for logs, varying in time since death, using negative exponential models and lag time models (Harmon et al. 1986(Harmon et al. , 2000)).The basic negative exponential model that includes only time since death (TSD) as a parameter had the form: where Y t is the density at time t, Y 0 is the initial density and k is the decomposition constant.The lag time model includes an additional parameter (n) to account for the time it takes for decomposer organisms to colonize and establish in a piece of newly created deadwood: where n represents the time lag in years.We incorporated TSF and diameter into Eqs. 1 and 2 to by altering initial wood density and/or decay rates (k) to produce a series of competitive decay models.We compared both linear and power-law functions to assess the effects of TSF and diameter on decay rate (k).We presented the linear function only when it was superior to the power-law function according to the AIC c .

Modeling of stand-level deadwood dynamics
We expressed the yearly input of deadwood biomass as the difference between the estimated biomass of living trees from the segmented regression models during the previous year and the current year.We estimated annual growth rates by measuring ring widths from growth years 2002-2007 of 20 living trees representative of the dominant cohort from the initial four stands sampled in 2008.Using these data, we then developed a diameter-dependent growth model based on models developed by Coomes et al. (2005) and G. Sainte-Marie (unpublished data).Our model included an additional parameter for time since fire (Eq. 3) where a was the effect of time since fire (TSF ), diam was the DBH and MaxRG, xo and xb were all constants (4.45, 9.30 and 0.77, respectively).Model constants were derived from Quebec provincial data by G. Sainte-Marie (unpublished data).One year of modeled annual growth was then applied to all stems in each plot and biomass of trees was recalculated.Both previous year biomass and current year biomass (previous year plus growth) were modeled using segmented regression.Current year live tree biomass was always less than previous years biomass plus growth, therefore the difference between the current year and previous year can be interpreted as tree mortality, or newly created deadwood biomass.Snag fall rate (SFR) in a given year was calculated as the difference of the current years snag biomass and the previous years snag biomass plus the inputs from living biomass (Eq.4) where m snags t is the mass of snags at time t and m inputs t is the mass of inputs from tree death.Note that the snag fall rate includes both snags transitioning to logs and trees that transition directly to logs.Decay constant of downed deadwood (k) was then calculated based on the negative exponential model of wood decay as the difference in the current year log biomass and the previous year's biomass plus inputs form snags minus logs buried in the organic layer (Eq.5) where m logs t was the mass of logs at time t, m snas.inputs was inputs from SFR and r burial t was the current year burial rate.We modeled three different scenarios for the rate of burial of logs on the assumption that burial rate is proportional to the rate of increase of soil organic layer.The yearly organic layer depth was estimated from the regression models and the yearly percent increase was determined as the difference be-tween current and previous years organic layer depth, divided by the previous years biomass, and multiplied by 100.The first scenario used the percent of biomass of logs buried each year equal to the percent increase of the organic layer.This scenario predicted much more buried log biomass than observed in the two measured stands.The second scenario used the percent increase of the organic layer divided by 5.This scenario predicted levels approximately between the two measured stands.The third scenario divided the increase by 10, which corresponded to the lowest value we observed in the two measure stands.Finally, we recalculated wood decomposition rates using a modified lag-time model, where exponential decay was permitted only after a lag period derived from the log decomposition time series.

RESULTS
General stand characteristics are presented in Table 1.Older stands had generally greater organic layer depth, lower basal area and, with the exception of the oldest stand, higher volumes of deadwood.In most stands nearly half of the logs were considered buried based on common designation [surface area of the log covered in .50%organic layer (Hagemann et al. 2009)].Additional sampling for buried deadwood using trenches in two plots [N50 (373 years since fire) and N6 (710 years since fire)] revealed an additional 110 m 3 /ha and 30 m 3 /ha of deadwood, respectively.
We found a significant relationship between height and DBH of trees using the Chapman Richards model (Table 2).The best model using the AIC c explained 60% of the variation of height in the data and did not include the effect of time since fire (TSF ) on any of the coefficients in the model.Both diameter and TSF affected living wood density.Live trees had a mean wood density of 0.40 g/cm 3 .Wood density of live trees increased with the log of TSF (GLM, b TSF ¼ 0.015, t ¼ 6.91, P , 0.001) and decreased with DBH (GLM, b DBH ¼ À0.0031, t ¼ À5.59, P , 0.001).

Changes of biomass with time since fire
We found strong support for breakpoints in the relation between TSF and biomass of live v www.esajournals.orgtrees, snags and logs (Davie's test p ¼ 0.007, p , 0.001, 0.007, respectively), and therefore used a segmented regression to model these variables.Live tree biomass declined at high rates (t ¼ À2.77, p ¼ 0.020; Fig. 3A) until 143 (SE ¼ 21) years after fire, after which biomass continued to decline but at a much slower rate (t ¼ À2.36, p ¼ 0.040).Initially, we did not observe a significant relationship between snag biomass and TSF.However, following the removal of sites N23 (a young site with a high number of snags) and Gau (an older site with a low number of snags), we found snag biomass significantly increased (t ¼ 2.79, p ¼ 0.027; Fig. 3B) until 255 years after fire (SE ¼ 43), at which point snag biomass significantly decreased (t ¼ À3.07, p ¼ 0.015).We chose to remove these two observations as they clearly did not fit the trend observed in the other sites and to allow us to fit a regression line that describes the majority of the sites for our stand dynamics model.We found a similar relationship between TSF and biomass of logs (Fig. 3C), although before the breakpoint, the increase in biomass was not significant (t ¼ 1.90, p ¼ 0.087).
However, following the breakpoint of 262 years (SE ¼ 64) since fire, the biomass of logs significantly decreased (t ¼ À2.26, p ¼ 0.047).Organic layer depth ( p ¼ 0.162) increased linearly with log (TSF; t ¼ 4.92, P , 0.001; Fig. 3D).We highlighted the trends of the site ''Pui'' in Fig. 3 as this site has many characteristics of a much older stand than our data indicated.For example, the relatively low biomass of live trees and the thick depth of the organic layer were much more consistent with an older stand or a stand that arose from a lower intensity fire.Unlike the other stands in the chronosequence, the age of this site was not verified by radiocarbon dating.Differences in stand characteristics of ''Pui'' may occur because this stand could be older or follow a different successional trajectory from a low intensity fire (Lecomte et al. 2006a).

Estimation of decay rate
We used a total of 488 cross-sections measured from 105 trees across 7 stands to estimate decay rate of logs.Decay models and AIC c values based on the negative exponential function are reported in Table 3 with parameter estimates in Table 4.The best model using the negative exponential (Table 3, model 8) indicated that as stands aged initial wood density increased (e.g., for every 100 year increase in stand age, initial wood density increased by about 1-3%) and decay rates decrease (e.g., for every 100 year increase in stand age, decay rates decreased by 2%).Decay models that included a lag before the onset of   6) had a lag of 6.6 years; for every 1 cm increase in diameter, initial wood density increased by 0.0022 g/cm 3 and for every 100 year increase in TSF, initial density increased by 0.0131 g/cm 3 and the decay rate (k ¼ 0.010) decreased from 14-6% with every 100 year increase of the mean TSF (x ¼ 242 years) over the first 1000 years.
Resampling of wood density of the same logs provided little evidence of decay in the initial two years.We were unable to observe a reduction in wood density in over half of the logs (10/19).In the remaining logs wood density was reduced by an average of 2.4% per year.Resampling of wood density of older logs (5-10 years since death) over a 2-year period revealed marginal evidence of wood decay.Two of the 18 logs could not be found for the re-measurement and again there was no detectable reduction in wood density of four of the logs.Statistics on the remaining 12 logs indicated a significant reduction in wood density (RM-ANOVA, F 1,13 ¼ 18.312, p , 0.001), with a reduction of an average 3.7% per year.However, limiting statistical v www.esajournals.organalysis to just those samples where mass loss was detectable may be biasing the data to predict a larger mass loss than actually occurring.
We also found that the best predictors of moss growth on re-sampled logs were the diameter of the log itself (LME, t 30 ¼ À3.41, p ¼ 0.002) and canopy openness (LME, t 30 ¼ 1.55, p ¼ 0.132).Although canopy openness was not significant, this model outperformed the model that included just log diameter (AIC c ¼ À299 vs. À292, respectively).

Transition of biomass from live trees, snags and logs
Tree growth rates (measured as ring widths) across the chronosequence decreased with increasing TSF and decreased with increasing diameter at breast height (Fig. 5A).Overall biomass of deadwood inputs decreased with TSF (Fig. 5B).The highest rates of deadwood input were in young stands during the time of highest individual tree growth rates and highest rates of decrease in live tree biomass.Snags fell faster in young stands (,143 years) compared to later in the chronosequence (Fig. 5C).
Under all three scenarios of burial rate, the predicted volume of buried deadwood fell within a range that could be expected based on the two trenches dug in two of the older stands (N50 and N6, shown in Fig. 6A as black circles).We used the negative exponential model to determine decomposition rates based on the inputs of log biomass from snags (Fig. 5C), the losses to the buried deadwood pool under the three above scenarios (Fig. 6A) and the observed biomass of logs present across the chronosequence (Fig. 3C).We found no difference in decomposition rates across the chronosequence under the three scenarios of burial rates (Fig. 6B).Our model predicted extremely high decomposition rates (k ¼ 0.544) during the senescence of the post-fire cohort early in succession, which then decreased to k ¼ 0.066 at 290 years post-disturbance, and then slowly increasing to k ¼ 0.113 as forests become more paludified.The average across the entire chronosequence was k ¼ 0.096.
Comparison of the two approaches (log vs. stand chronosequences) and the two decay models (negative exponential vs. negative expo- v www.esajournals.orgnential with a time lag) revealed that all four predict the highest rates of decomposition in the youngest stands, this rate decreases with stand age with the log chronosequence (Fig. 7).However, the stand chronosequence approach predicted slightly higher rates of decomposition in the oldest stands (.1000 years) than the intermediate stand ages (100-300 years).We also found that the two approaches which incorporated a time lag have more similar decomposition rates than the approaches that did not include a time lag (Fig. 7B).

Carbon storage
Aboveground carbon stored in live coarse woody biomass (.5 cm diameter) ranged from 16 Mg/ha to 54 Mg/ha (Fig. 8A).The younger stand (,100 years) had on average twice the levels of stored C than the older stands.The carbon stored in dead coarse woody biomass ranged from 2.3 to 10.3 Mg/ha with older stands having much higher amounts of stored C than younger stands in the dead pool.When estimates of buried coarse woody biomass were included (Fig. 8B), carbon storage in woody materials in black spruce remained relatively constant, ranging between 40 and 45 Mg/ha, throughout the chronosequence.

DISCUSSION
In this study, we demonstrated that reduced stand productivity due to paludification reduced the rates of deadwood input.Compared to mature forested stands, older paludified stands have less living woody biomass where a smaller proportion of that biomass dies each year resulting in less deadwood created each year.Conversely, we observed the highest rates of Table 5. AIC c and estimates of coefficients of models that include a lag time.v www.esajournals.orgdeadwood input early in succession during the senescence of the post-fire cohort.Our three approaches to estimate wood decay (resampling individual logs, comparing logs with different time since death and stand-level modeling) corroborated a delay in the onset of decomposition in individual logs and the highest decay rates when deadwood input is the highest.We also demonstrated that the highest rates of organic layer growth occurred when the canopy breaks up, which led to high rates of log burial.

Model number
As the buried wood accumulated in the organic layer it greatly increased the amount of stored woody carbon in these forests.
The cause of retrogression in these black spruce ecosystems is largely attributed to the paludification of the forests (Simard et al. 2007).Paludification is characterized by the reduction of large tree density.We found that in black spruce this reduction, measured as all living woody biomass with a diameter of greater than 5 cm, was much higher during the senescence of the post-fire cohort than later in succession even after using a log-transformed stand age.These results are similar to the models of Harper et al. (2005), who demonstrated that live tree density stabilized between 95 and 164 years after fire depending on site productivity.Other authors have characterized the relationship using a constant linear regression model comparing tree biomass to a log-transformed stand age (Lecomte andBergeron 2005, Simard et al. 2007).Both approaches predict high rates of tree loss during the senescence of the post-fire cohort; however, the non-constant regression parameter describes a relatively higher rate of loss during this earlier in succession and slower rate 140-160 years after fire.Regardless of the particular model to describe this decrease in tree density, the result is ;50 Mg/ha of deadwood biomass being created over an approximately 50-year period.
When predicting the quantity of deadwood being created each year based on changes in live tree biomass, it is important to include estimates of annual growth.When we included estimates of tree growth into our calculations of total deadwood input, we found that young stands had the highest growth rates, adding 2.1 Mg/ha of biomass per year, or 105 Mg/ha of biomass over this 50-year high deadwood input.We believe our estimates of growth are conservative, as they were based solely on the oldest trees.Our estimates were slightly less than the estimates of Simard et al. (2007), from the same forest type Fig. 4. Proportion of mass remaining with increasing time since death for five ages of forests using the best (lowest AIC c from Tables 3 and 5) negative exponential model (solid lines) and the time lag model (dashed line).
v www.esajournals.org(Fig. 9).Differences in our growth estimates and those reported in Simard's study could be explained by the fact that they included smaller trees which usually have larger annual radial growth rates.Therefore, even greater quantities of deadwood may have been created than predicted by our models.This could result in a greater biomass of stored carbon in the form of buried deadwood.
Our models predicted high inputs of deadwood biomass (155 Mg/ha) during the senescence of the post-fire cohort; however, these amounts do not readily reflect the observed aboveground biomass of deadwood.We did not observe maximal levels of aboveground deadwood biomass (22 Mg/ha) until 100 years following the senescence of the post-fire cohort when rates of deadwood input were lower.We believe that this extreme input of deadwood biomass rapidly leaves the stand through decomposition or is buried by the growing organic layer and incorporated into belowground carbon stocks.

Decay rates (k)
The majority of deadwood is lost from the forest through wood decay.We were able to improve our estimates of decay rates by incorporating time lags into exponential models of decay.In these forests, measurable decomposition began approximately 6-7 years following tree death.Many studies have found that the single negative exponential curve (Olson 1963) does not adequately describe wood decomposition (Grove et al. 2000, Harmon et al. 2000, Yatskov et al. 2003, Strukelj et al. 2013).Trees are thought to decompose slowly at first because of high initial moisture content, decay resistant heart-wood, or time required for decomposer organisms to become established (Grier 1978, Harmon et al. 2000, Yatskov et al. 2003).
However, not all other studies on black spruce have observed this lag in decomposition.For example, Boulanger et al. (2011) reported decay v www.esajournals.orgrates of fire-killed black spruce (k ¼ 0.027-0.036)similar to ours (k ¼ 0.037) using the simple exponential model which did not include a lag period.Although little evidence of decay is reported in black spruce snags (Angers et al. 2012), we believe that when a tree dies and remains standing as a snag for several years that the time lag associated with decomposition of logs is negated.We attribute the lack of a decompositional lag to preconditioning the snag by colonization of decomposers and loss of the high initial moisture content following death.Similarly, fire could also negate this decompositional lag as moisture content quickly decreases Fig. 6. (A) cumulative buried deadwood using three rates of burial relative to organic layer (OL) growth rates (Points represent additional buried deadwood found in trenches dug in two of the older sites) and (B) Decomposition rate (k) from negative exponential models under the same three burial rates.Fig. 7. Comparison of the time-series and stand-level modeling approaches for estimating decay rates using (A) negative exponential models, and (B) lag-time models.Both models use a burial rate proportional to a fifth of the organic layer growth rate.v www.esajournals.org in the fire.Therefore, lag-time decay models are limited to deadwood that comes into contact with the ground shortly after tree death.

Effects of diameter and stand age on decay rates
Wood decay was also affected by log diameter when using the negative exponential models.Diameter effects on wood decomposition rates have a strong theoretical basis where low surface-to-volume ratios are thought to prolong the colonization process and reduce sapwood-toheartwood ratios (Foster and Lang 1982, Harmon et al. 1986, Rayner and Boddy 1988).Diameter effects have been frequently reported on relatively small diameter pieces (Abbott and Crossley 1982, Barber and Van Lear 1984, Erickson et al. 1985) or between large and small size classes (Brown et al. 1996).However, many studies did not report significant diameter effects on woody decay rates (Foster and Lang 1982, Marra and Edmonds 1996, Bond-Lamberty et al. 2002, Grove et al. 2009).There are only a few studies which demonstrate this relationship across larger diameter pieces.For example, He ´rault et al. ( 2010) found a large effect of diameter on decay rates in Amazonian forests of French Guiana.A study in the Russian boreal by Tarasov (1999) [cited in Yatskov et al. (2003)] found a negative correlation between size and decay rates.Chambers et al. (2000) found a significant effect of DBH on decay rates, although this relationship explained little of the variation in the data (r 2 adj ¼ 0.10).v www.esajournals.orgHowever, we did not detect large effects of diameter on decay rates using the lag time model.Instead, time since fire (stand age) had the biggest influence on decay rates.In the range of stand ages (60-710 years) studied, the decay rate decreased 64% from k ¼ 0.18 in young stands to k ¼ 0.11 in older stands.This reduction in decay rates corresponds to low organic layer decomposition rates in paludified stands, which have been explained by colder and wetter environments, partly through the dominance of Sphagnum spp. on the forest floor (Gordon andShugart 1989, Lavoie et al. 2005).These wetter environments in paludified stands could lead to higher wood moisture levels, which is also known to reduce wood decomposition rates (Progar et al. 2000, Barker 2008).Bond-Lamberty et al. (2002) found that wetter black spruce sites had lower decomposition rates when wood moisture was above 43%.

Burial of deadwood
During forest paludification extremely slow decomposition rates of sphagnum mosses results in a rapidly accumulating soil organic layer (Fenton et al. 2010).We expected that the highest burial rate would correspond to the highest rates of accumulation of the soil organic layer, typical of the younger stands in our chronosequence.As expected, we observed effects of log diameter on burial rates; however, the rate of moss accumulation on logs was in older open stands, after the period of high organic layer growth.The rapid colonization of logs by mosses in these older, more open stands can be explained by shifts in moss composition.In the older stands (N50 and N6), the moss community on the logs was dominated by both Sphagnum spp.and humicolous species, whereas in the younger stands (N23, N18) only humicolous species were present on the logs.Furthermore in the oldest stand examined (718 years since fire), 4 different species of Sphagnum spp.were present whereas in the next youngest stand (373 years since fire) only a single Sphagnum sp. was present.The increasing dominance of Sphagnum mosses in older stands is also responsible for a rise in the water table (Lavoie et al. 2005), resulting in cool wet conditions effectively halting wood decomposition (Hagemann et al. 2009).Ultimately, burial of logs by the organic layer must related to the growth of the organic layer itself, even if we see higher rates of growth on individual logs in older forests when rates of organic layer growth have slowed.

Stand-level modeling of deadwood dynamics
We developed a simple deadwood model for black spruce forests by modeling changes in live and dead woody biomass over time providing an overall picture of deadwood dynamics throughout succession.This model allowed us to predict annual inputs of deadwood biomass and independently verify decay rates calculated from the results of decay models based on a time-series of logs.Our model included several assumptions of aspects of deadwood dynamics we were unable to directly measure.The model predicted extremely large inputs of deadwood during the senescence of the post-fire cohort; far exceeding the relatively low levels of dead biomass observed in the stands during this same period.Although this period was also predicted to also have the highest wood decay rates in both the time-series analysis and the deadwood model, the actual decay coefficient predicted by the deadwood model was much higher (k ¼ 0.54 compared to k ¼ 0.18).The average decay rate predicted across the entire chronosequence was k ¼ 0.096).These decay rates are much higher than any previously reported, especially in this period during the senescence of the post-fire cohort.For example, (Boulanger and Sirois 2006) reported decay rates for fire-killed black spruce of k ¼ 0.02, and the highest rates Brais et al. (2006) reported for all boreal species was in trembling aspen with a decay rate of k ¼ 0.06.These rates are, however, similar to ours (k ¼ 0.037) when using similar time-series approaches and considering just the negative exponential models of decay.Possible discrepancies between the stand-level model and time-series approaches could come from imprecision in growth rate or burial rates.For example, if we did not consider annual growth at all, decay rates would be reduced by a third.However, the trees are growing from year to year, and our estimates of growth are similar to other studies from the region (Simard et al. 2007).The other possibility is that a larger proportion of deadwood is consumed by the growing organic layer.This also seems unlikely, as increasing the burial rate by a factor of 5 or 10 had little effect on decay rates.Furthermore, the organic layer is not so deep that huge mass of deadwood would not be detectable by surface transect.
The final possibility is that this period during the senescence of the post-fire cohort does indeed experience exceedingly high wood decay rates.The large and constant input of deadwood could be supporting a unique decay community capable of consuming the deadwood created.Clearly, this period of high deadwood input and decay represents a unique period of deadwood dynamics across all forested ecosystems.These forests are losing an estimated 150 Mg/ha of deadwood biomass over 50 years.As we do not observe anywhere close to this quantity of deadwood in the forest, it must be decaying rapidly resulting in exceptionally short residence times on the forest floor.This will have direct impacts on deadwood-associated organisms, carbon budgets and forest management decisions.

Effects of deadwood processes on the carbon budget and forest management
Forest management strategies in these forests have focused on either clearcutting to restart the successional process or partial cutting to advance forest succession (Bergeron et al. 2002).To restart succession, clearcutting must be accompanied by severe disturbance of the forest floor by either mechanical scarification or prescribed burning (Lafleur et al. 2011) or paludification of the forest will continue.Alternatively, partial cutting has the greatest potential for maintaining carbon stock consistent with unmanaged forests if partial cutting is done during the maximal biomass phase before the senescence of the post-fire cohort (Fenton et al. 2013).We observed rapid rates of wood decomposition and provided evidence of high rates of burial during this period, thus residence time of surface woody biomass is at its shortest.Partial cutting during this period would also result in reduced woody biomass that would otherwise be rapidly decomposing.Therefore, reducing the amount of deadwood during or just prior to the senescence of the post-fire cohort should have relatively short-term effects on deadwood associated organisms as deadwood naturally leaves the system rapidly during this period compared to later in succession when deadwood has much longer residence times in the forest.Later in succession when the forests have entered a retrogressive state, deadwood and carbon stocks are at a more stable equilibrium and anthropogenic disturbance at this point would have longer lasting effect on the natural dynamics of deadwood and the organisms and processes that depend on deadwood.
Our study supports the idea that large pools of deadwood exist under the surface of the organic layer.These pools of buried deadwood have been shown to decay at dramatically slower rates and are thought to represent a longer-term pool of carbon than deadwood found on the surface (Moroni et al. 2010).We also demonstrated a possible mechanism of this process in black spruce stands undergoing paludification.When we combine our buried deadwood model with other measures of woody biomass, we found that these paludified stands maintained a relatively constant biomass of woody carbon following the senescence of the post-fire cohort.This highlights the importance of belowground deadwood biomass in forest carbon models (Kurz et al. 2009).Although forests that enter this retrogressive state continually decease in productivity and lose woody carbon aboveground to encroaching sphagnum mosses, they are able to continually store more woody carbon belowground as deadwood dynamics remains relatively constant for thousands of years.Long-term carbon storage in the form of buried wood will thus be promoted by leaving old forests intact.Harvesting of old forests, on the other hand, has the potential to reset succession and return the stand to a productive state.Management of these old forests then becomes a trade-off between goals of long-term carbon storage and wood production.

Fig. 1 .
Fig. 1.Map of the samplings sites in northwestern Quebec.Grey area on inset map denotes the Clay Belt region which is prone to paludification.

Fig. 2 .
Fig. 2. Conceptual model of deadwood dynamics in black spruce forests and the framework of the deadwood dynamics model.Numbers in brackets associated with arrows correspond to equations presented in the methods section.

Fig. 3 .
Fig. 3. Relationship between time since fire (years) and (A) live woody biomass, (B) snag biomass, (C) log biomass, and D) organic layer depth.Lines represented the fitted lines of segmented regression (A, B and C).All lines are significant ( p , 0.05) except for logs (B) before the breakpoint ( p ¼ 0.087).Grey boxes represent 95% confidence interval of the break-points.Black points in graph B represent points removed from analysis.The xaxis is on a log scale in graphs A-C.Points labeled ''Pui'' are from a site that appears older than date determined by oldest tree.

Fig. 5 .
Fig. 5. Results from deadwood model based on (continuation of Fig. 5 caption) concept in Fig. 2, (A) modeled growth across tree diameters and TSF, (B) inputs into the deadwood pool, and (C) proportion of snags that falls in a given year across the chronosequence.

Fig. 8 .
Fig. 8. Carbon storage of coarse woody biomass: (A) showing the relative levels of woody carbon stored in live trees, snags and logs with diameters above 5 cm, and (B) the same figure with the inclusion of buried logs using the most conservative scenario of burial rate from Fig. 6A.

Fig. 9 .
Fig. 9. Comparison of growth estimates across the chronosequence from our study and a more detailed analysis on forest productivity by Simard et al. (2007).Data does not include mortality between years.

Table 1 .
General stand characteristics of the 15 stands along the chronosequence.

Table 2 .
Chapman-Richards model parameter estimates for tree height based on DBH.The Chapman-Richards formula is yðtÞ ¼ að1 À e ðb 3 DBHÞ Þ c :

Table 4 .
Estimates of coefficients for wood decay models in Table3.Values in boldface represent the best model (i.e., lowest AIC c ).

Table 3 .
Wood decay models and AICc based on the negative exponential function.

Table 6 .
Estimates of coefficients for wood decay models in Table 5.Values in boldface represent the best model (i.e., lowest AIC c ).