© INRA, EDP Sciences, 2010

1. Introduction

There is an increasing interest to better estimate the amount of carbon sequestrated by forests. The contribution of the aboveground component can be assessed from tree volume and wood density as well as their variations in space and time. There is widespread evidence that forest productivity has increased in the last few decades in Europe (Spiecker et al., 1996) as a consequence of global climatic change (Kahle et al., 2008). Studies on carbon sequestration usually consider specific wood density to be constant (De Vries et al., 2006), and thus, there is little information about its possible long-term variations (Briffa et al., 1998; Bergès et al., 2000).

Mean ring density is known to vary according to (i) the genetic provenance (Hannrup et al., 2004; Rozenberg et al., 2001); (ii) tree age and growth intensity (Bouriaud et al., 2004; Jyske et al., 2008) (the effect of cambial age –ring number counted from the pith– is strong and conspicuous among tree species, while that of ring width depends on wood structure, Zobel and Van Buijtenen, 1989); and (iii) environmental conditions such as local climate (Briffa et al., 2002) and site fertility conditions (St-Germain and Krause, 2008). This has two implications: first, historical variations in mean density are plausible in the context of global climatic change; second, detecting these variations requires controlling all other sources of variation.

Figure 1 Geographical positions of the 13 sites sampled.

Investigations reporting long-term trends with growth-rate and cambial-age effects filtered out are insufficient for determining historical variations. An exception is provided by Bergès et al. (2000) for the case of sessile oak (Quercus petraea Liebl.) in a temperate climate. Whereas a long-term radial growth increase was reported, a decrease in wood density, assessed after allowing for changes in ring width and cambial age, was also evidenced. Given that the ring structure in sessile oak shows light earlywood and denser latewood, larger rings should imply higher ring density. The control of ring width thus revealed a downshift in the growth-density relationship. However, Bergès et al. assessed the abovementioned changes at a low temporal resolution (four values for one century) and without accounting for site variations. Decreases in maximum latewood density, corrected for the effect of cambial age only, were also found in conifers growing at a medium elevation in North America (Conkey, 1988) and all over the boreal treeline (Briffa et al., 1998). Despite methodological differences, the reported decreases in wood density show some common features, and this may raise some concern about both forest carbon sequestration and wood properties.

Norway spruce is an important coniferous species in the European forest resource. Therefore, it is essential that an analysis of its historical wood density variations is carried out. The objective of this study was to determine whether the mean ring density of Norway spruce growing in temperate conditions has followed the aforementioned negative trend. Spruce wood density is known to decrease with ring width (Olesen, 1976). As there has been an increase in the radial growth of Norway spruce in north-eastern France (Badeau et al., 1996), a decrease in wood density was expected. The specific objectives of this study were therefore (i) to account for variations in mean ring density that could be attributed to different growth rates, cambial ages, sites, and trees; and (ii) to check whether there exists a residual historical trend for wood density after accounting for the effects on mean ring density introduced by different growth rates, cambial ages, sites and trees, and to depict its chronology over the last century. In order to accurately separate the individual effects of cambial age and calendar date on the mean ring density, we propose an original method based on the paired-plots sampling design (Bontemps et al., 2009) and a combined statistical modelling approach.

2. Materials and methods

2.1. Sampling

The paired-plots sampling method (Bontemps et al., 2009) was used in order to detect a long-term change in the mean ring density of Norway spruce. Stands with sufficiently different ages were selected in each site according to this sampling scheme. This allowed us to compare the density of rings of identical cambial ages formed during different calendar years in the same site fertility conditions, and therefore to discriminate the effects of cambial age and date.

We sampled 13 sites in the Lorraine region of North-eastern France (France, Fig. 1), as this area represents an important part of the French forest resource in Norway spruce (Picea abies Karst. L.) (Seynave et al., 2005).

Sites were located in even-aged stands planted in State forests at various altitudinal positions (230 to 1 000 m a.s.l.). Different fertility classes assessed from soil organic layers and topographic situations were sampled.

While two stands per site were initially targeted, a majority of sites offered the opportunity to sample stands of three distinct ages, and thus to further reduce the correlation between age and calendar year. In total, 35 stands were sampled from thirteen sites. The analysis was restricted to dominant trees, as they correspond to the potential trees of the final harvest and are of primary interest regarding wood properties. Following Duplat and Tran-Ha (1997), three dominant trees were sampled per stand. In total, 105 trees were sampled. Tree ages across the 35 stands showed a bimodal distribution, with a separation of around 80 y between the younger and older trees. For these two classes, the mean values for ring width, mean ring density, and cambial age are given in Table I and Figures 2 and 3.

2.2. Wood density measurements

For each tree, we sampled one disc at breast height. From each disc, four radial strips were sawn for X-ray wood density measurements. According to the methodology described by Polge and Nicholls (1972), the strips were sawn at right angles to the fibre direction following the delimited radii. The strips had a constant thickness of 2 mm, and were conditioned to standard moisture content (ca. 12%). They were exposed for 4 h to X-ray radiation using Kodak Industrex AA 400 film and standard electrical conditions (accelerating tension = 7.5 kV; flux intensity = 12 mA), with the source 2 m away from the samples. The resulting radiographs of the samples and the plastic calibration wedge were digitized and analysed using specific image analysis software (Mothe et al., 1998, CERD software) providing ring widths once ring limits are semi-automatically located. X-ray density values must be calibrated with the gravimetric wood density (Freyburger et al., 2009). The weight of each sample was measured to the nearest 0.001 g and its volume was calculated from the thickness and area of the samples.

Annual ring density and ring width series were cross-dated by visual superimposition. In the following, mean ring density refers to mean ring density measured by the CERD software.

Average ring width and ring density were computed for each cambial age over the four strips.

The first 7 rings of each tree (juvenile wood, Olesen, 1976) were excluded so as to remove the strong negative cambial age effect on mean ring density close to the pith. Altogether, 5 212 rings were processed in this study.

2.3. Data analysis

The objective was to estimate the effect of calendar year on mean ring density series, after allowing for the effects of different growth rates, cambial ages, sites, and trees. Mixed-effects models efficiently handle between-individual variations in longitudinal data and have successfully been used for wood density analysis (Guilley et al., 1999; Le Moguedec et al., 2002). Both site (Tab. I) and tree differences (not shown) were detected in the data and included in the model. Site and tree differences were considered as random variations around parameter average estimates and were assumed to be Gaussian (Lindström and Bates, 1990).

We assumed a linear structure for the mean ring density model. Hence, the total variation of mean ring density could be split into a fixed predictive part, a random part associated with the site and tree levels of the sampling design, and a residual part (Hervé, 1999).

Parameters were estimated by maximum likelihood (ML) using the lme function from the nlme package (Pinheiro et al., 2009) of the R software (version 2.9.2). Selection of the best model was based on an analysis of the residuals, the Akaike information criterion (AIC), and the likelihood ratio Chi-2 test (LRT) for nested models (Pinheiro et al., 2009). All variables were corrected to zero mean in order to adjust parameters more easily; with the exception of the intercept, this operation did not change fixed parameters and random standard-deviation estimates.

The first step of this modelling phase was to place emphasis on the fixed part of the model. As depicted in Figure 3, the relationships between mean ring density to cambial age and ring width were both strong and curvilinear.

For each effect, the highest correlations with mean ring density were obtained using the square root of cambial age (r = 0.48) and a hyperbolic form of the ring width (r =  −0.6), which were introduced into the model (Eq. (1)).

Figure 2 Mean ring density distributions for (a) the whole dataset, (c) older trees, and (e) younger trees, and ring width distributions for (b) the whole dataset, (d) older trees, and (f) younger trees.

Figure 3 Mean ring density plotted against (a) cambial age and (b) ring width.

Model fitting was initiated by adding random effects of tree and site to the intercept alone. In order to avoid unnecessary duplication, the names and meanings of model variables for indices are given after equation (4). Equation (1) is as follows: $\begin{array}{c}\\ \mathit{W}{\mathit{D}}_{\mathit{ijk}}& \mathrm{=}& \mathit{m}\mathrm{+}\mathit{a}\mathrm{·}\sqrt{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}\mathrm{+}\mathit{b}\mathrm{·}\frac{\mathrm{1}}{\mathrm{1}\mathrm{+}\sqrt{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}}& \mathrm{Fixedpart}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{j}}& \mathrm{Site}\\ & & & \mathrm{Random}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{ij}}^{\mathrm{\prime }}& \mathrm{Tree}\\ & \mathrm{+}& {\mathit{\epsilon }}_{\mathit{ijk}}& \mathrm{Residuals}\mathit{.}\end{array}$(1) An analysis of the residuals suggested a remaining effect of cambial age, which was taken into account by an additional linear effect of the variable as shown in equation (2): $\begin{array}{c}\\ \mathit{W}{\mathit{D}}_{\mathit{ijk}}& \mathrm{=}& \mathit{m}\mathrm{+}\mathit{a}\mathrm{·}\mathit{C}{\mathit{A}}_{\mathit{ijk}}\mathrm{+}\mathit{b}\mathrm{·}\sqrt{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}\mathrm{+}\mathit{c}\mathrm{·}\frac{\mathrm{1}}{\mathrm{1}\mathrm{+}\sqrt{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{j}}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{ij}}^{\mathrm{\prime }}\\ & \mathrm{+}& {\mathit{\epsilon }}_{\mathit{ijk}}\mathit{.}\end{array}$(2) An analysis of the relationship between ring width and mean ring density for different cambial ages revealed that the relationship became weaker as the cambial age increased, showing that an interaction between cambial age and ring width was necessary (Eq. (3)): $\begin{array}{c}\\ \mathit{W}{\mathit{D}}_{\mathit{ijk}}& \mathrm{=}& \mathit{m}\mathrm{+}\mathit{a}\mathrm{·}\mathit{C}{\mathit{A}}_{\mathit{ijk}}\mathrm{+}\mathit{b}\mathrm{·}\sqrt{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}\mathrm{+}\mathit{c}\mathrm{·}\frac{\mathrm{1}}{\mathrm{1}\mathrm{+}\sqrt{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}}\mathrm{+}\mathit{d}\mathrm{·}\frac{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{j}}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{ij}}^{\mathrm{\prime }}\\ & \mathrm{+}& {\mathit{\epsilon }}_{\mathit{ijk}}\mathit{.}\end{array}$(3)Once the fixed part was assessed, random-effect parameterisations of variables in addition to pure between-site and between-tree variations (accounting for intercept variations in the model) were investigated. Given that our final objective was to detect an effect of calendar year, between-tree variation in the effect of age was not tested as it was likely to capture the historical trend. However, between-site variation was tested (cambial age assumed to be the same across cohorts in a single site). Variations in the effects due to ring width were investigated in order to detect between-site and between-tree variations. No random effect was tested on the interaction between cambial age and ring width. Following this procedure, the site and tree effects were formulated as expressed in equation (4). The resulting model (M1) is expressed as follows: $\begin{array}{c}\\ \mathit{W}{\mathit{D}}_{\mathit{ijk}}& \mathrm{=}& \mathit{m}\mathrm{+}\mathit{a}\mathrm{·}\mathit{C}{\mathit{A}}_{\mathit{ijk}}\mathrm{+}\mathit{b}\mathrm{·}\sqrt{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}\mathrm{+}\mathit{c}\mathrm{·}\frac{\mathrm{1}}{\mathrm{1}\mathrm{+}\sqrt{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}}\mathrm{+}\mathit{d}\mathrm{·}\frac{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{j}}\mathrm{+}{\mathit{\beta }}_{\mathit{j}}\mathrm{·}\sqrt{\mathit{C}{\mathit{A}}_{\mathit{ijk}}}\\ & \mathrm{+}& {\mathit{\mu }}_{\mathit{ij}}^{\mathrm{\prime }}\mathrm{+}{\mathit{\gamma }}_{\mathit{ij}}^{\mathrm{\prime }}\mathrm{·}\frac{\mathrm{1}}{\mathrm{1}\mathrm{+}\sqrt{\mathit{R}{\mathit{W}}_{\mathit{ijk}}}}\\ & \mathrm{+}& {\mathit{\epsilon }}_{\mathit{ijk}}\end{array}$(4)where WD_ijk is the mean ring density of ring k in tree i at site j, CA_ijk is the cambial age, RW_ijk the ring width, (m, a, b, c, d) are the parameters of the fixed part, (μ_j, β_j) are the random parameters associated with the between-site variations, (, are the random parameters associated with the between-tree variations, and ε_ijk is the residual term assumed to be Gaussian and independent of the random-effects (Lindström and Bates, 1990).

The last modelling step focused on the detection of a historical trend conditional upon the effects of growth rate and cambial age. The general form of the model is therefore defined by equation (5): $\mathit{W}{\mathit{D}}_{\mathit{ijk}}\mathrm{=}\mathit{f}\mathrm{(}{\mathit{C}{\mathit{A}}_{\mathit{ijk}}\mathit{,R}{\mathit{W}}_{\mathit{ijk}}}^{\mathrm{)}}\mathrm{+}\mathit{g}\mathrm{(}{\mathit{t}}^{\mathrm{)}}$(5)where f(CA_ijk, RW_ijk) is the model expressed in equation (4) and g(t) is the effect of calendar year. This effect was described with three different time resolutions (centennial, decadal, annual), which are respectively defined as follows.

• (i)

  A linear function of date for detecting a long-term trend (model M2) $\mathit{g}\mathrm{(}{\mathit{t}}^{\mathrm{)}}\mathrm{=}\mathit{p}\mathrm{·}\mathit{t}$(6)where p is the date parameter and t is the expression of the date modified in order to set 1900 as a reference date.

• (ii)

  A cubic spline function of the date for detecting decadal variations along the time period covered (Bontemps et al., 2009). This function is a continuous piecewise polynomial. Different node intervals were tested (10, 20, and 30 years). Using AIC comparisons, a significantly better fit was obtained with a 20-years internode (model M3): $\begin{array}{ccc}\mathit{g}\mathrm{(}{\mathit{t}}^{\mathrm{)}}& \mathrm{=}& {\mathit{d}}_{\mathrm{1}}\mathrm{·}\mathit{t}\mathrm{+}{\mathit{d}}_{\mathrm{2}}\mathrm{·}{\mathit{t}}^{\mathrm{2}}\mathrm{+}{\mathit{d}}_{\mathrm{3}}\mathrm{·}{\mathit{t}}^{\mathrm{3}}\\ & & \mathrm{+}\sum _{\mathit{i}\mathrm{=}\mathrm{1}}^{\mathrm{4}}{\mathit{p}}_{\mathit{i}}\mathrm{·}{\left[\begin{array}{}\max \end{array}\mathrm{(}{\mathrm{0}\mathit{,t}\mathrm{-}\mathrm{20}\mathrm{·}\mathit{i}}^{\mathrm{)}}\right]}^{\mathrm{3}}\mathrm{+}{\mathit{p}}_{\mathit{m}}\mathrm{·}{\left[\begin{array}{}\min \end{array}\mathrm{(}{\mathrm{0}\mathit{,t}}^{\mathrm{)}}\right]}^{\mathrm{3}}\end{array}$(7)where t is the expression of the date modified in order to set 1900 as a reference date, and (d₁, d₂, d₃, p_i, p_m) are the spline parameters to be estimated.

• (iii)

  An annual effect of the date treated as a factor (model M4), so that one average variation of the mean ring density per year was estimated. This allowed us to characterize year-to-year variations.

For each model, we computed the decomposition of variation using the methodology developed in Hervé (1999).

In order to quantify the importance of the mean ring density trend detected (conditional upon the effects of variability in ring width, cambial age, tree, and site), the detected trend was compared to the raw ring density variation over the sample and the variation that would result solely from growth rate change, as estimated from the model. We first selected all rings with a cambial age below 40 years in order to have comparable mean cambial ages for older and young trees, and then calculated differences in mean ring density and mean ring width between young and old trees. The calculated difference in mean ring width represented the shift in radial growth over the century and was used to predict (Eq. (4)) the partial wood density difference attributable to changes in growth rate, all other things being equal. We therefore obtained three values of wood density change: a raw ring density change, a ring density change attributed to the shift in ring width, and the historical shift remaining when controlling for all effects.

3. Results

3.1. Relationships of mean ring density to cambial age and ring width

Overall, the average ring density was 439.5 kg m^-3 with a standard deviation of 62.7 kg m^-3, i.e. a coefficient of variation of 14% (Fig. 2). The corresponding figures for ring width were 3.3 mm, 1.6 mm, and 48%, respectively. The breakdown is given in Table I.

3.2. Modelling the effects of age and ring width on mean ring density

Table II shows the main features of model M1. As expected, all parameters except the intercept were found to be significant since the variables were centred around a zero mean. Random variations indicated the variability of parameters according to site and tree. For the particular case of the intercept, the standard deviation associated with tree was twice as high as that associated with site. This indicated that the between-tree variability was higher than between-site variability in the sample.

Figure 4 Residuals of model M1 plotted against (a) cambial age and (b) ring width.

As shown in Figure 4, the distribution of residuals was satisfactory.

The effects of ring width and cambial age were thus considered accurately represented. A hyperbolic function of ring width and its associated positive parameter were consistent with the decrease in the mean ring density with the ring width. Parameters associated with the effect of cambial age showed a sharp increase in mean ring density near the pith, which then became more gradual.

The fixed part of model M1 explained 33.5% of the variation in mean ring density (Tab. III). The total random variation due to the tree and the site represented 43.9% of the total variation. This random part was mainly located at the tree level (31.7%) and was linked to the variation in the intercept. The remaining random variation came from variation in the intercept according to site, while the between-tree and between-site variations of cambial age and ring width parameters remained weak. The residual variation remained high (22.6%), and an analysis of variance (data not shown) revealed that the year-to-year variation observed in the residuals was significant (p-value  <  0.0001).

3.3. Historical changes in mean ring density

Results for models M2 to M4 including an effect of calendar year are shown in Table II. Regardless of model, all parameters of the fixed part were stable, with a slight difference between models M3 and M4 in the parameter associated with cambial age. Moreover, the standard errors for the parameters of the fixed part and the standard deviations for both random parts were also constant.

The year-to-year effect of date was highly significant (1132 point difference in AIC) as compared to the linear trend (LRT p-value  <  0.0001). Decomposition of the variation also showed that the contribution of the fixed part of the model (Tab. III) increased from M1 to M4. The fixed part represented around 33% of the variation for M1 and nearly 45% for M4, consistent with variations in the AIC. There was no change in the importance of the random effect of site. In contrast, the variability associated with the total random effect of tree decreased by 5% from M1 to M4. The inclusion of a factorial function of year in the fixed part of the model strongly diminished the residual variability.

Historical trends in the mean ring density at different time resolutions are displayed in Figure 5. First, there was an average decrease in the mean ring density of up to 18 kg m^-3, corresponding to 5% of the mean ring density. However, this date effect was weak (p-value = 0.0348). The polynomial form showed that the effect of calendar year was somewhat constant in the first part of the twentieth century. In the second part, there was a steep decline in the mean ring density of more than 30 kg m^-3 with an increase at the end of the century. The factorial effect of date followed the polynomial form: there was no visible effect until 1960, after which there was a steep decrease in mean ring density. Moreover, there was very high year-to-year variability. From one year to the next, the mean ring density varied by up to 60 kg m^-3, reflecting the high standard deviation of the mean ring density in our data. From the different forms of the date effect introduced in models M2 to M4, we concluded that during the twentieth century, the mean ring density was not constant but decreased by 5% during this period.

A direct comparison of the average ring density between young and old trees revealed a decrease of 52 kg m^-3, which was higher than the evolution found in the previous model M2. As the average ring width increased by 1.68 mm (from 2.65 mm for old trees to 4.33 mm for young trees), the partial prediction of model M1 based on this difference corresponded to a decrease in the mean ring density of 32 kg m^-3. The trend estimated with model M2 displayed a diminution of 18 kg m^-3, all other things being equal. The reduction in the mean ring density observed when comparing young and old trees therefore had two origins: an increase in the ring width and an additional even reduction with ring width and cambial age held constant.

Figure 5 Long-term trend, decadal fluctuations, and inter-annual variations in mean ring density throughout the twentieth century.

4. Discussion

4.1. Determinants of mean ring density in Norway spruce

The relationship between ring width and density was strong and negative, as has already been reported in previous studies on Norway spruce (Mäkinen et al., 2007; Molteberg and Hoibo, 2007). Earlywood has a larger lumen and thinner cell walls than latewood (Decoux et al., 2004), an increase in ring width leads to a smaller proportion of latewood (Mäkinen et al., 2002) and accounts for the observed density relationship.

We also found a positive relationship between mean ring density and cambial age. This is consistent with previous studies on Norway spruce (Cao et al., 2008; Pape, 1999). For most tree species, there is a general decrease in the ring width with cambial age in mature wood (Gutierrez Oliva et al., 2006; Pape, 1999). Given the negative relationship between ring width and wood density, the observed effect of cambial age may thus be a consequence of the decrease in ring width with cambial age. However, in accordance with our results, models expressing wood density as a function of both ring width and cambial age in Picea abies show a significantly positive and direct effect of cambial age (Ikonen et al., 2008; Wilhelmsson et al., 2002). One interpretation of this effect is based on the increasing proportion of latewood with cambial age (Mäkinen et al., 2007). In addition, an increase in cambial age leads to greater cell wall thickness and a larger cross-sectional area of the cell (Lundgren, 2004a).

In addition to these effects, between-tree and between-site variations were also evident (Tab. III). Variation among individual trees is well documented (Zobel and Buijtenen, 1989), and the present range of variation for mean ring density is consistent with this observation: it accounted for 20 to 25% of the total variation in wood density across our sample. This is usually interpreted as an effect of genetic variation (Gerendiain et al. 2007; Hannrup et al. 2004; Rozenberg et al., 2001). However, information on genetic origin is lacking in our sample. Site variation amounted to 10% of the total variation (Tab. III), and also affected the relationship between cambial age and mean ring density, in agreement with a study by Raiskila et al. (2006).

4.2. Historical trends in mean ring density

The relationships of mean ring density to cambial age and ring width are usually assumed to be constant with time. However, we hypothesized that these relations may have been modified in the changing global climate. These modifications may thus result in a historical trend, a prediction which was tested in the modelling phase.

We found a moderate but significant decrease of the mean ring density throughout the twentieth century, totalling 18 kg m^-3 (5%). This range is consistent with that found in a previous study on oak (Bergès et al., 2000).

Conkey (1988) and Briffa et al. (1998) assessed the trend in conifer species. However, they did not take the effect of ring width on wood density into account, which could have led to a bias in their results. By controlling for the growth rate, we found a less pronounced decline closer to the estimate reported by Bergès et al. (2000). We hypothesize that this decrease may be quite general since for trees it has been assessed for diverse wood structures and in different growth conditions around the world.

A 70% increase in radial growth was found when comparing young and old trees (Tab. IV), consistent with a previous study on Norway spruce in the Vosges mountains (Badeau et al., 1996). According to the relation between mean ring density and ring width (Mäkinen et al., 2007), this corresponded to a decrease of 32 kg m^-3. We conclude that the mean ring density has decreased throughout the twentieth century because of two phenomena: an effect of the increase in ring width (and changes in ring textures) and an historical trend which was highlighted by keeping all other sources of variation constant. An increase in radial growth coupled with a decrease in mean ring density was also found in oak in the same proportion by Bergès et al. (2000). However, as the structure of oak wood is different from Norway spruce, the total effect on mean ring density in this previous study is a combined consequence of the radial growth increase and a pure effect of the calendar year. Moreover, in our study, the historical variation in mean ring density is similar in magnitude to the effect of site variations across the sample (assessed from the random site effect in Tab. II).

Our modelling approach allowed us to extract evidence of decadal and annual variations in mean ring density. Three chronological periods can be identified: (i) from 1900 to 1950, when ring density decreased slowly; (ii) from 1950 to 1980, when a steep decrease in the ring density occurred reaching a ring density reduction of 30 kg m^-3 relative to 1900; and (iii) from 1980 to 2000, when a final increase occurred from –30 kg in 1980 up to –15 kg m^-3 reached in 2000. This pattern is remarkably similar for different species (Bergès et al., 2000), locations, and environments (Briffa et al., 1998; Conkey, 1988), and supports the hypothesis of a widespread phenomenon in trees. However, a major difference in our study is that we were able to produce a chronology and to estimate the date effect up to the end of the century after allowing for the differences in the cambial age, ring width, site, and trees.

Bergès et al. (2000) suggested that the decrease in wood density could have been a manifestation of the heartwood formation process. The border between heartwood and sapwood in Norway spruce is not characterized by a natural change in colour (Longuetaud et al., 2006) and the content of extractives in Norway spruce is low. It could be regarded as negligible contributors of wood density and the extractives removal could thus be neglected (Jaakkola et al. 2005). Differences in wood density occur mainly for fresh wood (Karenlampi and Riekkinen, 2003) and wood density does not vary at low moisture content (Taylor et al., 2002). Moreover, chronologies show that the trend is neither linear nor monotonous, and is thus in conflict with the heartwood formation process.

The chronology at an annual resolution (Fig. 5) was characterized by high year-to-year variability (standard deviation of around 20 kg m^-3). Annual variations in ring density are known to be mainly explained by variations in environmental factors such as temperature (Fritts, 1976).

4.3. The potential role of environmental factors

The effect of an increase in carbon dioxide on wood structure is poorly documented, and the available results appear to be contradictory in Norway spruce. Hattenschwiler et al. (1996) reported an 8% increase in earlywood density under elevated CO₂ concentrations for spruce seedlings. Similar increases of 4.8 and 2% in earlywood and latewood densities, respectively, were also obtained by Kostiainen et al. (2004) by doubling the CO₂ concentration under natural growing conditions in a mature plantation. However, Kostiainen et al. (2009) also found a 9% decrease in the latewood density by doubling the CO₂ concentration, together with anatomical variations in the whole ring – a decrease in cell wall thickness and increase in tracheid diameter. The latter result is consistent with a report by Kilpelainen et al. (2005) on Scots pine. There was no apparent relationship between the reported mean ring density chronology and the centennial increase in CO₂ (IPCC, 2007).

The reported effects of precipitation on wood density have varied in sign among studies (Gagen et al., 2006; Jyske et al., 2010; Yasue et al., 2000). Moreover, precipitation did not vary significantly during the last century (Moisselin et al., 2002). The role of precipitation is thus unclear.

Dendrochronological studies usually show a strong positive effect of temperature on wood density in Norway spruce (Wang et al., 2002; Wimmer and Grabner, 2000). The reported 1950–1980s decrease in ring density thus appears to be in conflict with the temperature stability during this period (Moisselin et al., 2002). Our observations are therefore consistent with the sensitivity loss of the mean ring density to temperature, which was pointed out by Conkey (1988) and confirmed by Briffa et al. (1998). However, a recent increase in mean ring density during a period of unusual warming (Moisselin et al., 2002) is also reported. To our knowledge, such an increase has not been found previously for mean ring density. As a particular result of the present study, we thus postulate a recent recovery in the temperature sensitivity of the mean ring density. This suggests a role for transient environmental influences during the period of divergence between wood density and temperature (Briffa et al., 1998; Conkey, 1988).

The deposition of nitrogen on the Earth’s surface has strongly increased in the second half of the twentieth century (Blake et al., 1999; Hastings et al., 2009), resulting in the eutrophication of most ecosystems (Vitousek et al., 2000). The effect of nitrogen availability on wood density has been reported as being negative based on fertilization experiments. In the study by Mäkinen et al. (2002), the observed decrease may be a consequence of the growth increase which was also reported. However, Lundgren (2004b) also found a decrease in wood density of 30 kg m^-3 while controlling for both ring width and cambial age. Lower wall thickness and higher cell diameter were observed in support of the observation (Lundgren, 2004a). Therefore, the nitrogen hypothesis may be compatible – in both the sign and timing – with the decrease in mean ring density reported here. It may also account for the previously reported loss in temperature sensitivity (Conkey, 1988).

5. Conclusion

The mean ring density of Norway spruce has declined during the twentieth century. This decrease has two origins: (i) an increase in radial growth and (ii) a further decrease in ring density after allowing for changes in the effects of ring width and cambial age as well as between-site and between-tree variations. This trend in the mean ring density is not linear, but shows both increases and decreases with a steep decrease in the second part of the century ending with a slight increase in the last 20 years. These trends in mean ring density show some consistency with the evolution of nitrogen deposition and are partially in accordance with the change of temperature. This is particularly noticeable in the last decades of the twentieth century and is consistent with the unprecedented global climatic change.

As our measurements permit analyses of within-ring density variations, the next step will be to separately analyse the wood density variations of earlywood and latewood. As latewood is known to be the most sensitive to environmental variations, the goal will be to investigate the relationships between environmental variables and latewood density throughout the last century.

Acknowledgments

The PhD grant of the first author was funded partly by the DGER (Direction Générale de l’Enseignement et de la Recherche) from the French Ministry of Agriculture and partly by the Lorraine region (France). We thank the French National Forest Service (ONF) for providing financial support through the Modelfor project. And we also gratefully thank the two anonymous reviewers for detailed and useful comments on an earlier version of the manuscript.

References

All Tables

All Figures

	Figure 1 Geographical positions of the 13 sites sampled.
In the text

	Figure 2 Mean ring density distributions for (a) the whole dataset, (c) older trees, and (e) younger trees, and ring width distributions for (b) the whole dataset, (d) older trees, and (f) younger trees.
In the text

	Figure 3 Mean ring density plotted against (a) cambial age and (b) ring width.
In the text

	Figure 4 Residuals of model M1 plotted against (a) cambial age and (b) ring width.
In the text

	Figure 5 Long-term trend, decadal fluctuations, and inter-annual variations in mean ring density throughout the twentieth century.
In the text