The use of soil nutrient balances in deriving forest biomass harvesting guidelines specific to region, tree species and soil type in the Netherlands

The substitution of biomass for fossil fuels in energy consumption is a measure to decrease the emissions of greenhouse gases and thereby mitigate global warming. During recent years, this has led to an increasing interest to use tree harvest residues as feedstock for bioenergy. An important concern related to the removal of harvesting residues is, however, the potential adverse effects on soil fertility caused by increased nutrient removal, relative to conventional stem-only harvesting. In the Netherlands this is a major concern, since most forests are located on poor sandy soils. To develop forest harvesting guidelines, we applied a mass balance approach comparing nutrient inputs by deposition and weathering with nutrient outputs by harvesting and leaching for various timber harvesting scenarios, including both stem-only harvesting and additional removal of tree tops and branches. A distinction was made in seven major tree species, six soil types (three sandy soils, loam, clay and peat soils) and nine regions, with clear variations in atmospheric deposition of phosphorus (P), calcium (Ca), magnesium (Mg) and potassium (K). For each region-tree-soil combination we calculated the maximum amounts that can be harvested such that the output of the nutrients Ca, Mg, K and P is balanced with the inputs. Results showed that at current harvesting rates, a negative balance of Ca, Mg, K or P is hardly calculated for the richer loamy to clayey soil types, while depletion can occur for the poorer sandy soils, particularly of P and K. Results are used to derive forest biomass harvesting guidelines, taking the uncertainties in the mass balance approach into account. The role of mitigating management approaches is also discussed.


Introduction
Forests are expected to be used more intensively in the next decades, with growing demands for timber and for use of biomass as a source of renewable energy (Mantau et al., 2010;Nabuurs et al., 2015). A key reason for the expected increase in the use of biomass, including wood biomass, for energy production in the coming decades is the need to reduce GHG emission, associated with using fossil fuels (e.g. European Commission, 2016). Biofuels based on forest biomass consist of harvesting residues, including tree tops and branches, which can be co-harvested with conventional timber (Evans and Finkral, 2009). While current forest biomass harvesting levels in Europe and the Netherlands are below annual increment (Schelhaas et al., 2014), a growing market for wood materials and biofuels will put increasing pressure on resource provisioning from forests (Mantau et al., 2010;Nabuurs et al., 2015). Whole tree harvesting is increasingly applied in Europe as the marketability of the woody residues has increased. This trend is in line with the policy of the EU and also of the Netherlands to reach a complete circular economy by 2050. By harvesting the entire above-ground part of the tree, the biomass removed by harvest can be 30% larger than compared to conventional stem-only harvesting (Pels, 2011, Aherne et al., 2012. The use of forest biomass as a source of renewable energy is highly debated, both in the scientific community and the public domain, especially in view of its sustainability in relation to climate and biodiversity (e.g. Searchinger et al., 2018). This debate focuses, however, on burning harvested wood for energy production, such as wood pellets that mainly come from tree stems of pulp wood quality or sawdust, otherwise used for wood products (Walker et al., 2015). It is argued that harvesting and burning wood biomass causes CO 2 emissions, just like any other carbon fuel (Hudiburg et al., 2019). The discussion of whole tree harvest versus stem wood only removal, however, focuses on the use of woody residues as biofuels, which otherwise would mainly be decomposed. Here, the key question refers to the sustainability of soil fertility, since an increased amount of nutrients is exported from forests (Raulund-Rasmussen et al., 2008, Iwald et al., 2013, Thiffault et al., 2014, Achat et al., 2015. Considering that element concentrations in needles, leaves and branches are much higher than in stem wood, whole tree removal may have a significant impact on soil quality as compared to stem only harvest (e.g. Aherne et al., 2012). Wall and Hytönen (2011) found no significant differences in the nutrient amounts in the forest floor and the foliar nutrient concentrations between whole-tree harvesting, with the needles left on the site, and a conventional stem-only harvesting treatment in a Norway spruce stand 30 years after clearcutting. However, based on a meta-analysis of published data worldwide, Achat et al. (2015) concluded that removing forest harvesting residues leads to an overall loss of soil fertility due to reductions in total and available soil nutrients, as a consequence of increased nutrient export, with a negative effect on the growth of subsequent forests.
Sustainable production of forests implies that the nutrient status should not decline over subsequent rotations, which depends on nutrient losses by forest export (harvest) and leaching on the one hand and inputs by deposition and weathering on the other hand. Apart from removal of nutrients through tree harvest, forests are most sensitive to nutrient losses through leaching in the developmental stages directly following harvest (Den Ouden et al., 2010). Such nutrient losses are particularly high for base cations in acidified sandy soils (De Vries et al., 1995). Increased harvest of biomass can thus reduce tree growth (Egnell, 2011;Kreutzweiser et al., 2008), lower nutrient concentrations in needles (Olsson et al., 2000), lower base cation pools in the soil (Olsson et al., 1996;Zetterberg et al., 2013) and increase soil acidity (Feller, 2005;Kreutzweiser et al., 2008;Achat et al., 2015). In addition removal of forest residues might reduce soil carbon sequestration (e.g. Clarke et al., 2015;Forsius et al., 2016). Based on a risk analysis of effects of whole-tree harvesting on site productivity, Wall (2012) concluded that phosphorus (P), and the base cations potassium (K), calcium (Ca) and magnesium (Mg) were priority nutrient indicators of site productivity, being reduced by enhanced nitrogen (N) deposition and soil acidification.
To guide, monitor and assess progress towards sustainable forest management (SFM), various criteria and indicators (C&I) have been developed (The Montreal Process, 1995). The Pan-European Forest Process on criteria and indicators for SFM (41 signatory countries), for example, developed six criteria including the maintenance and appropriate enhancement of (i) forest contribution to global carbon cycles, (ii) forest ecosystems health and vitality, (iii) productive functions of forests, (iv) biological diversity in forest ecosystems, (v) protective functions of forests (notably soil and water) and (vi) other socio-economic functions and conditions. These six criteria are made more concrete with 27 related quantitative and 84 qualitative indicators. Unsustainable use of nutrients may affect nearly all of those functions, but there are no criteria specifically related to it. Hence, governments and non-governmental entities have developed a range of governance tools, such as policies, regulations, certification schemes and guidelines, to ensure that forest biomass harvesting is sustainable (Stupak et al., 2007;Stupak et al. 2013). Most guidelines are, however, based on qualitative assessments.
In this paper, we present a quantitative nutrient balance approach, comparing inputs and outputs for different combinations of tree species, soil types and regions in the Netherlands, being the basis of forest biomass harvesting guidelines in view of sustainable nutrient use (see also Akselsson et al., 2007a;Forsius et al., 2016). The nutrients include P, Ca, Mg and K, of which the availability is limited due N deposition and soil acidification. Sustainability in terms of soil conditions was defined as a situation in which no long-term reduction in P, K Ca and/or Mg pool is allowed (Aherne et al., 2012;Akselsson et al., 2007aAkselsson et al., , 2007b). The idea is that mining of these nutrients implies a decrease in the concentrations of total N and P and of exchangeable base cations (Ca, Mg and K) in the mineral topsoil, being key soil fertility indicators in view of forest growth (Thiffault et al., 2014;Forsius et al., 2016). Actually Thiffault et al. (2014) also included total N in the mineral topsoil as a key indicator but there is no risk for N depletion in the Netherlands due to high N deposition caused by intensive animal W. de Vries, et al. Forest Ecology and Management 479 (2021) 118591 husbandry. The focus is on sandy soils, since these soils are particularly sensitive to acidification and loss of base cations (e.g. De Vries et al., 1995, De Vries et al., 2014, 2017a, and the major part of Dutch forests is located on those soils.

Nutrient balance approach
Nutrient fluxes that were included to assess the balance of soil nutrient inputs by atmospheric deposition and mineral weathering and outputs by forest harvesting and leaching (Fig. 1). The supply of nutrients through seepage streams or through flooding of forests with river or stream water was neglected since the majority of Dutch forest is located on well-drained sandy soils were seepage and flooding does not occur and the balance was thus calculated according to: Nutrient balance = Deposition + Weathering -Leaching-Harvest removal (1) The critical element removal (in kg ha −1 yr −1 ) by forest harvesting was set equal to the net input of P, Ca, Mg and K, implying no loss of available soil nutrient, according to: Critical nutrient removal = Deposition + Weathering -Leaching (2)) The critical wood removal (in m 3 ha −1 yr −1 ) was calculated by dividing the critical element removal for the most restrictive element (i.e. lowest critical element removal) by the element concentration in the harvested timber (in kg m −3 ), divided in stems and branches. Based on the deposition, weathering and leaching, we thus calculated how much wood (stems or stems plus branches in case of whole tree removal) can be harvested for different combinations of seven tree species, six soils and nine regions. The seven tree species were birch, beech, oak, Douglas fir, Norway spruce, Scots pine and larch. The six soils were three sandy soils (coarse poor sand, moderate poor sand and mineral rich or loamy fine sand), loess, clay and peat soils. The relationship between the different soil clusters and soil types, as distinguished on the 1: 50,000 soil map of the Netherlands, is given in De Vries et al. (2019). The nine distinguished geographical regions, differing in soils and/or typical habitats, are given in Figure S1. The focus was on sandy soils in the regions 2-5. because most Dutch forests occur on those soils and because sandy soils are most sensitive for depletion of Ca, Mg and K.

Calculation method
The removal of nutrients through wood harvest was calculated with the GrowUp model (Bonten et al., 2016). GrowUp is a tree growth model that calculates the uptake, retention and removal of nutrients per tree species per year based on the amount of biomass and the nutrient contents in stems, branches, roots, leaves and needles (See Figure S2 for a schematic representation of the model). With GrowUp, the net removal of Ca, K, Mg and P was determined for each combination of tree type, growth class (growth forecast) under two harvesting scenarios: (SOH) stem only harvest and (WTH) whole tree harvest at final felling, which included removal of needles for coniferous trees (Scots pine, Douglas fir and Norway spruce), combined with stem wood harvest at thinning.
The biomass increase calculations (in m 3 ha −1 yr −1 ) were based on growth curves (per tree species, growth class and harvest and thinning regime) that in turn were based on yield tables (Jansen et al., 1996). The mass of stem wood was calculated by multiplying the volumetric growth (in m 3 ha −1 yr −1 ) with the stem density (kg m −3 ). The mass of the branches, leaves or needles and roots was calculated on the basis of biomass expansion factors, originating from Vilén et al. (2005). The harvest of wood (ratio of thinning, final felling) is also based on those yield tables. For this study, only the ratio of branches to stem wood and the mass of needles at final felling was used, since other biomass components (especially roots) were assumed to remain in the forest and thus not included in the calculations. GrowUp also contains default values for nutrient concentrations in the biomass compartments based on literature research, but these values were replaced by measurements that were carried out for each of the considered tree species.
The nutrient removal by forest harvesting was calculated as the annual average removal over the rotation period of the relevant tree species, as used by Jansen et al. (1996), with a maximum of 100 years. The rotation periods given by Jansen et al. (1996) are 150 years for beech, 120 years for oak and pine, 100 years for Douglas, 90 years for birch and 80 years for larch and Norway spruce. The nutrient concentrations were assumed to be constant over the rotation period since average harvest levels over a whole rotation were calculated. These average nutrient concentrations were calculated by assessing the total nutrient and wood exports through thinning and final felling with varying proportions of wood compartments (bark, sapwood, heartwood) over time within a rotation (De Vries et al., 2019).

Nutrient levels in stems and branches
Nutrient levels in stems and branches were based on measurements in each of the seven tree species. Details on the sampled locations, the sampling and the processing of nutrient data in bark and wood (heartwood and sapwood) from stems, coarse branches and fine branches are given in De Vries et al. (2019). The calculated nutrient levels in stem wood and branch wood are given in Fig. 2. The same nutrient levels were used for all soil types, as samples from poor to rich sandy soils did not indicate clear differences in nutrient concentrations between the three categories of sandy soils.
For both stems and branches, N concentrations are highest, followed by Ca and K, while the P concentrations are relatively low. There are clear differences between tree species. For example, the levels in stems in oak and beech are higher than those in the other tree species. For Norway spruce, some levels (Ca, Mg) are higher than in the other coniferous trees and birch. In the branches, the calcium levels are highest for Norway spruce, followed by oak.

Input of nutrients by deposition
Calcium, potassium and magnesium: The spatial variation in the total (sum of wet and dry) atmospheric deposition of the base cations Ca, K and Mg was based on combined measurements (wet deposition) and model calculations ( Figure S3).
Phosphorus: The total deposition of phosphorus (P) is less accurate than for cations for the Netherlands, because dry deposition data are missing. Since 1992, however, P concentrations have been measured at 16 measuring stations of the RIVM rainwater monitoring network. Water quantities were also measured from which the wet deposition was derived by multiplying the precipitation amount by the P concentration. The wet P deposition of all measuring stations over the period 1992-2015 showed no trend and was on average around 0.08 kg P ha −1 yr −1 . Literature data show an average ratio near 2.3 for throughfall P-deposition / wet P-deposition with a range of 1.5-3.0 (Table 1). Throughfall might, however, overestimate the total P-deposition due to P leaching from forest crowns (Parker, 1983) and therefore a somewhat conservative average value of 2.0, with a range of 1.5-2.5, seems reasonable. An average value of 2.0 is in accordance with global estimates by Mahowald et al. (2008) and Vet et al. (2014). This ratio has been used as an average, which amounts to an average total P deposition of 0.16 kg P ha −1 yr −1 . Based on the variation between measuring points, a total P-deposition was calculated per 5 km × 5 km (Fig. 3), which was subsequently averaged for the nine geographical regions ( Figure S3).

Supply of nutrients by weathering
Weathering rates for the six distinguished soils (poor sand, moderately rich sand, rich sand, loess, clay and peat), have been derived from weathering experiments for sandy soils (De Vries, 1994) and for loess, clay and peat soils ( Van der Salm et al., 1998). Weathering of cations depends on the amount of weatherable minerals in the soil and increases with a decrease in pH. In allocating weathering rates, a generic reference pH was used per soil type. The weathering rates for a depth of 1 m of soil were based on (i) a revised evaluation of the data in Van der Salm et al. (1998), resulting in an update for loess-and clay soils and (ii) literature research into weathering rates based on field research and weathering models (Hyman et al., 1998;Klaminder et al., 2011;Yang et al., 2013;Starr et al., 2014;Johnson et al., 2015). Weathering rates for peat soils were based on the clay content of these soils. The results used are given in Table 2.
As with deposition, it also applies to weathering that relatively little is known about the supply of available P by weathering. Newman (1995) provides a range of 0.04-0.2 kg P ha −1 yr −1 for Europe based on a literature study. In our calculations, we used an average value of 0.1 kg P ha −1 yr −1 .

Leaching of nutrients
The leaching of the base cations Ca, K and Mg, is mainly linked to the leaching of sulphate (SO 4 ) and nitrate (NO 3 ) originating from acid deposition (cf. De Vries and Leeters, 2001). Leaching is the residual term from the mass balance and equal to deposition plus weathering minus forest harvesting plus a decrease in the readily available exchangeable base cation stock (see Fig. 1). The latter applies as long as there is an exchangeable supply of base cations that (partially) neutralizes the acid input. In fact, this is hardly or no longer the case in Dutch sandy soils. A re-sampling of Dutch sandy soils under oak in 2015 showed that in all forest stands the exchangeable base saturation has declined from around 10-20% in 1990 to around 5-10% in 2015 (De Vries et al., 2017a). In other words, there is currently hardly or no decrease in exchangeable base stock in sandy soils anymore and leaching is equal to deposition plus weathering minus uptake. Since we developed harvesting guidelines in view of uptake, we cannot include a value for uptake. We thus calculated leaching as a fraction of the deposition plus weathering, occurring mainly outside the growing period. This fraction was set at 0.5 (50%). Using this fraction caused a leaching rate of Ca, K and Mg that is quite comparable to a leaching rate that was estimated by multiplying measured average Ca, K and Mg concentrations levels in shallow groundwater (see Table 3) with an estimated annual average precipitation surplus (described below). Table 3 gives average measured Ca, K and Mg concentrations in shallow groundwater (at approximately 1 m below ground surface) at 18 locations in 2017-2018 and predicted annual average concentrations of Ca, K and Mg in water leaving the soil system. The predicted levels were calculated by dividing the calculated leaching of Ca, K and Mg, estimated at 50% of their deposition and weathering, by a 30-year annual average precipitation surplus for major forest types. Measurements were generally 5-15% higher than predictions, indicating that the leaching might have been slightly underestimated.
The 30-year annual average precipitation surplus was estimated as precipitation (30 year average) minus interception minus evapotranspiration, being the sum of soil evaporation and tree transpiration. The interception was set at a fraction of the precipitation, i.e. 0.40 for dark coniferous forests (Norway spruce and Douglas fir), 0.30 for light coniferous forests (Scots Pine) and 0.25 for deciduous forests (birch, beech, oak and larch), based on Rijtema and De Vries (1994). The evapotranspiration level was set at 350 mm yr −1 independent of tree species, based on the upper value of a range of 250-350 mm yr −1 , given by Roberts (1983). The upper value was used considering an elevated evapotranspiration in response to an increased precipitation in the last decades. Using a 30 year average precipitation of 850 mm, being the mean value for the period 1980-2010, we thus used a longterm average precipitation surplus for forests on sandy soils of 290 mm yr −1 for deciduous forest (birch, beech, oak) and larch, 240 mm yr −1 for light coniferous forest (Scots pine) and 160 mm yr −1 for dark pine forest (Douglas fir and Norway spruce). We neglected the variation in runoff between regions since the impact was considered insignificant when considering the uncertainty in interception evaporation of the different tree species.
Leaching of P (kg.ha −1 yr −1 ) was calculated by multiplying the above mentioned precipitation surpluses (m 3 .ha −1 yr −1 ) with an estimate of the annual average concentration of P in water leaving the soil system (mg P .l −1 ), because P is strongly buffered. We used an average P concentration of 0.02 mg of P l −1 based on the average value in shallow groundwater at 18 locations, as mentioned above.

Nutrient balances
Nutrient balances were calculated for 7 (tree species) × 6 (soil types) × 9 (regions) combinations to derive specific forest biomass  harvesting guidelines for each combination. The focus was sandy soils, because most Dutch forests occur on those soils and because sandy soils are most sensitive for depletion of Ca, Mg and K. The calculated Ca, K, Mg and P balances are illustrated for all tree species on clay, loess and peat soils (Fig. 4) and on sandy soils, split into coarse poor sand, moderate poor sand and rich fine sand (Fig. 5). As average growth during the rotation period in m 3 ha −1 yr −1 we used 4.5 for birch, 6.0 for oak, 8.0 for beech and Scots pine, 10.0 for larch, 11.0 for Norway spruce and 12.0 for Douglas fir. Results are also given for all sandy soils in the four sand regions (Fig. 6) and for each of the seven tree species (Fig. 7). In each figure, nutrient balances are given for stem wood only harvest (SOH) and whole tree harvest (WTH) at final felling. The figures show the inputs (deposition and weathering) on the positive part of the y-axis and the removal (harvest and leaching) on the negative part of the y-axis. The resulting balance term can be both positive and negative. A negative balance is on the positive part of the y-axis (red), being the extra supply of Ca, K, Mg or P that is required to maintain the Ca, K, Mg or P contents in the soil. The reverse applies to a positive balance, which is on the negative part of the y-axis (green), being the Ca, K, Mg or P surplus that increases the Ca, K, Mg or P contents in the soil. Averaged over all tree species, there was no negative balance for Ca or Mg, neither for SOH nor for WTH for the clay, loess and peat soils, but a negative balance (i.e. depletion of soil stock) for K on clay and peat soils and for P on all soils. For P the net decrease in the soil P pool was much larger for WTH than for SOH (Fig. 4). For the three sandy soils there was no negative balance for Mg, neither for SOH nor for WTH, while this was always the case for P when considering the average result of all tree species. As with clay, loess and peat soils, the negative P-balance, i.e. depletion of soil P, was much higher for WTH than for SOH. For K and Ca no depletion was calculated under SOH, but it occurred for K on all soils and for Ca for the poor sandy soils when applying WTH (Fig. 5). A similar result was calculated for the sandy soil regions. It shows that in all regions, there was on average no depletion of Mg, both for SOH and WTH, while there was always depletion for P, except for SOH in region 5. For K there was also depletion in all regions for WTH and for Ca, depletion occurred in region 2 and 5 for WTH. When removing only stem wood, no Ca or K deficits occurred in any region for SOH (Fig. 6). Fig. 7 shows a negative Ca balance for beech, oak and Norway spruce, both for SOH and WTH. For K, a negative balance was calculated for beech and oak for both SOH and WTH, while WTH also shows a negative balance for Douglas fir and Norway spruce. A negative P balance was calculated for beech and oak for SOH, and for all tree   (Fig. 7). Figure S4 shows two more specific examples, namely the balances for Scots pine on poor sand and oak on rich sand in the four sand regions. A noticeable difference is that pine on poor sand had hardly ever a negative balance for Ca, K and Mg, while this was the case for Ca and K for oak even on rich sand. This has to do with the higher demand from oak, which is only partly compensated by a higher weathering rate   (4) and South (5). for stem only harvest and whole tree removal. W. de Vries, et al. Forest Ecology and Management 479 (2021) 118591 from a richer soil. For both sites, we calculated a negative P-balance, being much higher for oak on rich sand than for pine on poor sand (in view of a higher demand), and positive Mg balances.

Use of nutrient balance results in forest biomass harvesting guidelines
The results of the balance calculations have been summarized in forest biomass harvesting guidelines for sandy soil by indicating the effects of harvest intensity for the seven studied tree species on poor, moderately rich and rich sandy soils, for each of the four sandy regions under SOH and WTH (De Vries et al., 2019). Since the results of the balance calculations did not show a strong decrease for loess, clay and peat soils and considering the limited area of forests on those soils, they were not included in the harvesting guidelines. Results in the Tables 4 and 5 refer to the Central sandy soil region while the Tables S1-S6 refer to the Northern, Eastern and Southern sandy soil regions (See Figure  S1). The tables show the nutrients of which the levels are estimated to decrease for a given combination of region, tree species and soil type.
The harvest levels indicated in the tables are the amounts of harvested wood averaged over the entire rotation, thus being the average growth level of the trees. Harvest levels are distinguished in low (4.5 m 3 ha −1 yr −1 ), average (7.5 m 3 ha −1 yr −1 ) and high (11.5 m 3 ha −1 yr −1 ), considering stem wood (diameter > 7 cm), using a 0.5 m 3 ha −1 yr −1 margin. The harvest levels that are achieved in practice in Dutch forests vary considerably and depend, among other things, on the increment and the functions/goals that the forest owner has on his forest and on the tree species. The increment of beech, birch, oak and Scots pine is on average clearly lower (mostly in a range of 4-8 m 3 ha −1 yr −1 ) than that of Douglas fir, Norway spruce or larch (mostly in a range of 8-14 m 3 ha −1 yr −1 ; see Fig. 8). An harvest level of 11 m 3 ha −1 yr −1 of deciduous trees and Scots pine is not considered applicable since the harvest intensity is larger than the actual growth. The average actual increment of larch, Norway spruce and Douglas fir is 9, 12 and nearly 14 m 3 ha −1 yr −1 , respectively. These three more productive tree species mostly have a production function, whereby most of the increment is ultimately harvested over the entire rotation. For those species, therefore, the harvest levels of 7-11 m 3 ha −1 yr −1 are particularly relevant.
Tables 4 and 5 show that for beech and oak, the Ca, K and P stocks are estimated to decrease at moderate and high (when applicable) growth levels, not only for WTH but also for SOH. At a low harvest level, there was no Ca depletion for beech and no Ca and K depletion for oak. In stands of Douglas fir and Norway spruce, soil nutrient stock can decrease at high growth level under SOH and also at a moderate growth level for WTH. In stands of birch, Scots pine and larch, the depletion

Table 4
Risks for depletion of nutrients (Ca, K and/or P) when harvesting stem wood only (SOH) in the Central sandy soil region of the Netherlands 1 .

Tree species
Poor sandy soils Moderate poor sandy soils

Rich sandy soils
Harvest intensity (m −3 ha −1 yr −1 ) Harvest intensity (m −3 ha −1 yr −1 ) Harvest intensity (m −3 ha −1 yr −1 ) 4.5 7.5 11.5 4.5 7.5 11.5 4.5 7.5 11.5 Elements are included when the input of nutrients (by deposition and weathering) is lower than their output (by harvest and leaching) at given harvest level, implying nutrient depletion. The abbreviation "na" means "not applicable" (irrelevant) since the harvest intensity is larger than the actual growth of the considered tree species. The parts in white imply that the inputs of all nutrients are either equal or more than the outputs.
was limited to P under SOH, except for K depletion in birch stands on poor sandy soils, while K depletion occurred in more stands under WTH (compare Tables 4 and 5). It is striking that a negative balance of Mg was not calculated in any of the combinations of tree species sandy soil types and regions (Tables 4 and 5 and Tables S1-S6). There were some differences between the regions, mainly due to differences in deposition, but they were not striking. Table 6 indicates in a qualitative sense the estimated uncertainties in the various input and output fluxes, and methods to reduce the uncertainties. The uncertainties are indicated in terms of low, moderate, considerable and high, in which estimated deviations from the reference, i.e. < 15% (low), 15-30% (moderate), 30-50% (considerable) and > 50% (high), value are based on an educated guess. The uncertainties in the balance terms are low for harvest removal, moderate to considerable for deposition and considerable to high for weathering and leaching, in particular for P.

Uncertainties in the calculated balances
The uncertainty in the element removal by harvest is on average low (< 15%) for all elements, since growth levels are well known and nutrient concentrations in stems, branches and leaves have been measured. Site-specific nutrient concentrations, however, often differ from the measured average concentrations. The uncertainty in the deposition is moderate for base cations (15-30%) and considerable for P-deposition (30-50%). This is because P deposition is measured at less locations than base cations, and because dry P-deposition is uncertain, which makes the total deposition estimate more uncertain. The uncertainty contribution of deposition is also dependent on the region. For example, it is relatively large for P and Ca in the southern region and relatively large for base cations in the coastal regions. The greatest uncertainty is the degree of leaching. This uncertainty is considerable for the bases (30-50%) and high for P (> 50%). In addition, leaching, in particular of base cations, makes a relatively large contribution to the total balance, comparable to the contribution of harvesting.
As with leaching, the uncertainty in weathering is considerable for Ca, Mg and K (30-50%) and high for P (> 50%). Regarding base cation weathering, quite a lot of research has been carried out, both in the Netherlands and internationally, and in general the results for sandy soils and loamy soils are fairly reliable. The assigned range of weathering rates of 250-600 mol c ha −1 yr −1 , based on weathering experiments for sandy soils, corresponds to literature data based on field research and weathering models as presented before. However, there is debate as to whether the current weathering rate is not much higher due to the lower pH values. Bergsma et al. (2016), for example, estimated the recent weathering rate for a sandy soil in the Veluwe region in central Netherlands based on mineral depletion in the topsoil over a period of 75 years. They found a loss of bases in that period of 95 kmol c ha −1 for a layer of 30 cm, which amounts to an annual weathering of about 1300 mol c ha −1 yr −1 . That is 2.5 times as high as the assigned weathering rate for rich sandy soils. However, the measured concentrations of bases in shallow groundwater give no reason to assume these high weathering rates. With regard to P weathering, there is hardly any literature available, which makes the estimate very uncertain. Newman (1995) provides a range of 0.04-0.2 kg P ha −1 yr −1 for Europe based on a literature study, implying an uncertainty of > 50% in the assigned value of 0.1 kg P ha −1 yr −1 .
An aspect that also causes uncertainty in the estimates, at least for the future, is the ongoing climate change, including an increase in temperature and precipitation, accompanied by an increase in CO 2 concentration. On one hand, this will likely increase the growth of Table 5 Risks for depletion of nutrients (Ca, K and/or P) when harvesting stem wood during thinning and whole tree harvest at final felling (WTH) in the Central sandy soil region of the Netherlands (see footnote for Table 4). forests (e.g. Kauppi et al., 2014;De Vries et al., 2017b), and thereby the removal of P, K, Ca and Mg from the soil, but it will also enhance the weathering rates of these nutrients (Qafoku, 2015). The increase in both temperature and precipitation will also affect the precipitation surplus, and thereby the leaching rate of nutrients, both by an increase and decrease, depending on the region. The overall effect of climate change on the nutrient balance is therefore hard to predict, but it increases the uncertainty in estimates.

Mitigating management approaches
There are various mitigating management approaches to reduce element decrease in the soil (De Vries et al., 2019), including nutrient addition and alternative harvesting methods. Adding nutrients is an effective method, but should be applied with care in view of potentially unwanted side effects. Adding P, if the P stock is insufficient, is a costeffective option without negative side effects since P is strongly bound to the soil. With respected Ca, K and Mg, it is best to use slowly dissolving fertilizers, such as rock powder, because liming is typically associated with a strong pH increase, enhanced decomposition and N release, leading to negative effects on biodiversity. However, adding rock powder is an expensive method as rock powder contains only a small proportion of the required nutrients per ton applied material. The costs would not be compensated by enhanced forest growth and harvest benefits. Alternatively, common agricultural Ca, Mg and K fertilisers could be used, but this should be done in low doses to avoid negative side effects.
Under whole tree harvesting, leaving branches on the soil for six months leads to a reasonable decrease in the harvest removal of K and P (up to 20%) and to a smaller extent of Mg and Ca (Staaf and Berg, 1982;Palviainen et al., 2004) due to leaching from the leaves. For Douglas fir, the measure is also fairly effective for Ca and Mg (relatively large needle mass compared to branch mass). Logistically, the measure is less attractive, because the harvesting of stems and branches takes place at different times. Leaving branch wood in the forest, implying stem wood only harvest, is particularly effective for Douglas fir and Norway spruce for P (around 40-60% less harvest removal), while the effects are relatively limited for birch, beech and oak (< 20% less harvest removal) and slightly more effective for Scots pine and larch (approx. 15-30% less harvest removal). The effect of longer forest rotation times and avoidance of clear-cutting while initiating rejuvenation quickly, may reduce the leaching after harvest but it is difficult to quantify these measures. When applying shelter wood cutting, as opposed to clear cutting additional nutrients, i.e. P, Ca, Mg and K, are captured by dry deposition, but this also holds for nitrogen that may acidify the soil and thereby enhance the leaching of base cations. Research on the fate of nutrients at forest locations where various mitigating measures are being implemented helps to better quantify the effectiveness of these measures.

Conclusions
Results of the nutrient balance calculations show that soils below beech and oak are most sensitive to loss of Ca, K and P stocks, even under stem only harvesting (SOH), except for Ca and K at low growth levels. Under whole tree harvesting (WTH), the K stock of beech and oak can also decrease at low growth levels. In stands of Douglas fir and Norway spruce, soil nutrient stocks can decrease at high growth level under SOH and also at a moderate growth level under WTH, while in stands of birch, Scots pine and larch, the depletion is mostly limited to P. A negative balance of Mg is never calculated. In order to improve the reliability of the forest biomass harvesting guidelines, additional data are needed, especially with respect to weathering and leaching of nutrients.

Declaration of Competing Interest
None.