Effects of stand features and soil enzyme activity on spontaneous pedunculate oak regeneration in Scots pine dominated stands – implication for forest management

A challenge in current forestry is adaptation of managed forests to climate change, which is likely to alter the main processes of forest dynamics, i.e. natural regeneration. Scots pine will probably lose some parts of its distribution area in Europe. However, two native oaks, pedunculate and sessile may maintain or expand the area of their occurrence in central Europe. The utilization of spontaneous (not initialized by foresters) oak regeneration in Scots pine stands for the creation of next generation stands is one of the adaptation methods to climate change. Many factors influencing pedunculate oak regeneration are well known, but there is a lack of knowledge on the relation between soil enzyme activity and the establishment and development of the species. The aim of the study was to identify the relationships among stand characteristics, herb species composition, soil enzyme activity and the establishment or recruitment of oak regeneration in Scots pine-dominated stands. The one of the most influential factors shaping the oak seedling count was dehydrogenase activity in the humus horizon. We found that plots without litter and fern cover had higher seedling density. The raspberry ground cover and birch crown projection area had a positive influence on oak seedling number. The factor indicating good conditions for high density of oak saplings was phosphatase activity in the organic horizon. The same enzyme activity but in humus horizon described conditions in which more numerous recruits were observed. The activity of soil enzymes can be used as the predictor of the establishment and advancement of oak regeneration but also could be seen as a new dimension of oak regeneration. The general density of spontaneous oak regeneration was not sufficient for the creation of new generation forest stands dominated by oak, but it is possible to use them as admixtures in new generation stands.


Introduction
Adaptation of managed forests to climate change is a challenge in current forestry (Lindner et al. 2014). Climate change is likely to alter the main processes of forest dynamics, i.e., tree growth, mortality, and regeneration (Seidl et al. 2016). Tree species will respond differently to climate change (Huang et al. 2017). Scots pine (Pinus sylvestris L.), which was promoted in previous centuries and covers large areas in the lowlands of northwestern and central Europe (Goris et al. 2007), will likely experience decreases in its distribution in Europe (Sáenz-Romero et al. 2017;Dyderski et al. 2018). However, other tree species may expand their ranges. Economically important tree species that may maintain or expand the area of their occurrence in central Europe include two native oaks, pedunculate (Quercus robur L.) and sessile (Q. petraea (Matt.) Liebl.) (Hanewinkel et al. 2013;Takolander et al. 2019). Both grow across the European lowlands and occur in many forest types (Eaton et al. 2016). Currently, Scots pine vitality has decreased in Poland (Report of Forest Status 2018) and in other countries. It has been observed in Europe that climate warming has made Scots pine stands more vulnerable to attack by bark beetles (Jaime et al. 2019) or mistletoe (van Halder et al. 2019). The changes in species distribution ranges could generate large-scale economic problems that are difficult to solve with predominant even-aged management because the bigger compositional adjustment could be made only at the end of production cycle and the cost of acceleration of such activities could be prohibitive (Schelhaas et al. 2015). One of the possible methods to increase the speed of response for problems connected to climate change is the utilization of oak spontaneous regeneration processes observed in current Scots pine stands for the creation of next generation stands or at least to enhance the species enrichment of current pine-dominated forest ecosystems of artificial origin (Galiano et al. 2010;Rigling et al. 2013).
Spontaneous regeneration of both oaks, silver birch (Betula pendula Roth), downy birch (Betula pubescens Ehrh.) or common beech (Fagus sylvatica L.) has been detected in many aging Scots pine forests (Zerbe 2002;Dobrowolska 2006;Kint et al. 2006;Gniot 2007;Goris et al. 2007). However, spontaneous natural regeneration established under the Scots pine layer is not often used, even though it can play a great role in the restoration processes of forest stands and enhances their adaptability to changing environmental conditions (Diaci et al. 2008;Vizoso-Arribe et al. 2014). The success of spontaneous oak regeneration depends heavily on animal activity. Numerous mammals and bird species contribute to oak dispersal by in European forests. Eurasian jays (Garrulus glandarius) and wood mouse (Apodemus sylvaticus) are the most important dispersers and hoarders of acorns in Europe. Most other species that forage on acorns are mainly seed predators or only occasional dispersers (Den Ouden et al. 2005). In the case of Scots pine dominated forest stands where older, seed producing oaks are sparse or absent (Mosandl and Kleinert 1998;Frost and Rydin 2000;Kurek and Dobrowolska 2016) jays seems to be a more important factor allowing colonisation of those stands by oaks than mice. Dispersal distances of acorns by jays are in the order of several hundred metres (Bossema 1979;Kollmann and Schill 1996) and much greater than dispersal distances by mice, which according to Den Ouden et al. (2005) ranged between 2.7 to 9.2 m with exceptional situations achieving 70 m. Eurasian jays prefer thick, deciduous forests but are also found in coniferous or mixed forests as well as in parks, gardens, and yards with plenty of mature trees (Vander Wall 2001). These birds hide acorns in hoards for winter time. The regeneration of oak establishes under a canopy of old pine and gradually grows up as more light becomes available in the gaps (Kint et al. 2004(Kint et al. , 2006. The regeneration of pedunculate oak has been intensively studied in Europe (Worrell and Nixon 1991;Paluch and Bartkowicz 2004;Harmer et al. 2005;Ligot et al. 2013;Annighöfer et al. 2015); however, knowledge of spontaneous oak regeneration in Scots pine-dominated stands is still limited. Moreover, most of the investigations have concentrated on the natural processes in aging pine stands (Kint et al. 2004). The results of these studies suggested that mixtures with oaks increase Scots pine growth and that those stands can be more productive (higher volume) than pure stands (Bielak et al. 2014;Steckel et al. 2019).
Many factors influencing pedunculate oak regeneration are well known (Kelly 2002;Annighöfer et al. 2015;Didenko and Polyakov 2018), but the knowledge on the impact of the soil enzyme activity on the establishment and development of the species is insufficient.
Enzyme activities are critical to ecosystem functioning, affecting nutrient transformation, carbon sequestration, and biogeochemical cycling of carbon, nitrogen, phosphorus and sulphur. Phosphatases catalyse the hydrolysis of ester bonds between phosphate and carbon compounds in organic substrates (Schneider et al. 2001). Urease activity may reflect the decomposition rate of nitrogen compounds in soil. The presence of asparginase in the substrate is responsible for the decomposition of organic nitrogen compounds, that supply nitrogen to plants. Dehydrogenases are responsible for oxidative reactions in soil (Wolińska and Stępniewska 2012). Enzyme activity correlates with soil physicochemical conditions, and nutrient availability (Tarafdar and Jungk 1987). Enzyme activities have been used as sensitive indicators of soil quality changes under the influence of management (Acosta-Martinez et al. 2014) and as an indicator of heavy metal contamination (Bojarczuk and Kieliszewska-Rokicka 2010;Yang et al. 2019).
Analysis of soil enzyme activity provides information on soil microbial status (Yang et al. 2019). One of the most important components of these organism are fungi that form mycorrhizae. Burke et al. (2011) showed that ectomycorrhizal fungal communities were positively correlated with most soil enzymes, including enzymes involved in C, N, and P cycling. Some observations linked the success of young tree regeneration to the presence of mycorrhiza, which could also be responsible for soil enzyme activity. Pedunculate oak seedling growth, budburst, survival, biomass and foliar nitrogen and phosphorus content depend on root colonization by ectomycorrhizal fungi (Newton and Pigott 1991;Egerton-Warburton and Allen 2001). Seedling growth is correlated with the number of mycorrhizal types (Newton 1991;García de Jalón et al. 2020). Moreover, oak regeneration surrounded by phylogenetically distant neighbours show increased abundance and enzymatic activity of ectomycorrhizal fungi in the litter (Yguel et al. 2014). The dependence between soil enzyme activity and various ecosystem characteristics described above suggests that they may interact directly or indirectly with oak regeneration establishment and advancement in our study area.
The aim of our study was to expand our knowledge on the spontaneous regeneration of pedunculate oak in Scots pine-dominated stands of different ages. In particular, we addressed the following questions: -Which stand characteristics influence the establishment and recruitment of oak spontaneous regeneration in Scots pine stands? -Does the presence of an herb layer and its species composition affect the density of oak regeneration? -Does soil enzyme activity influence the establishment and advancement of new oak generation? -Which factors influence the density of oak recruits belonging to different quality classes?

Study area
The study was conducted in northeastern Poland in the Masurian Lake District. Although it is the lowland part of Poland, the landscape is very diverse with many hills and lakes due to the influence of repeated glacial events, especially by the Baltic glaciations. In this region, the continental climate clashes with the Atlantic climate; thus, typical temperate zone forests grow in the Masurian Lake District. The average temperature in July is 18°C, the average temperature in January is − 4°C, the annual precipitation ranges from 500 to 634 mm, and the growing season lasts 190-200 days (Lorenc 2005).
The study was carried out in 2013-2015 in managed forests located in the Nowe Ramuki Forest District (138 m a.s.l.; 20°30′ E, 53°47′ N). The Nowe Ramuki Forest District is part of the great complex of the Napiwodzko-Ramucka Forest. The forest cover of the study area is high (67.9%), much higher than the average for the Masurian Lake District (Plan of Forest Survey in Nowe Ramuki Forest District). The main tree species is Scots pine, which occupies 91% of the forest area. Other important trees are pedunculate oak (4%) and silver birch (3%). The big share of Arenosols soils may suggest that in the primeval forest in this region the share of broadleaved trees (Quercus, Betula, Carpinus) was substantial. Most forests are old (average stand age is 83 years). Pine stands older than 100 years represent 41% of the area of the forest district. Clear-cutting (areas less than 4 ha) and planting are the usual methods used for Scots pine forest stand regeneration. However, spontaneous oak regeneration (originating from animal acorn dispersal) is frequently observed 1-2 years after clear-cutting, even if there were no oak trees in or nearby the former stand (foresters observations from the Nowe Ramuki Forest District).
The investigation was conducted in Scots pine stands, and the total area of our study covered 90 ha. The age of the Scots pine stands ranged from 26 to 140 years. All investigated stands were regenerated by planting on two site types: poorer (mixed coniferous forest site type, where the plant community Querco roboris-Pinetum is commonly observed) and richer (mixed deciduous forest site type where patches of Tilio-Carpinetum calamagrostietosum are abundant). These site types refer to the following soils: Dystric Albic Brunic Arenosols and Dystric Brunic Arenosols (WRB 2015).

Data collection
To investigate the influence of stand characteristics on oak regeneration, we collected data from 13 pinedominated stands of different ages (Additional file 1: Table A2). There are numerous pine stands in the study region that grow under very similar conditions to our studied stands, but have no or very sparse oak regeneration. The most likely factors responsible for this difference are limitations in oak dispersal related to the availability of seed sources and the activity of seed dispersing animals. Investigation of these factors would require a different experimental design and should be conducted at a larger spatial scale, and of course would require identification of potential seed sources. Our study focuses on other questions, namely the role of local stand conditions and soil enzyme activity in oak regeneration establishment. We chose a pine stand with abundant oak regeneration under the assumption that seed rain density there was not as limiting a factor as other factors considered in our experiment. This assumption might immediately raise the question of whether seed rain density limitations were comparable in all these stands and whether other factors not included in our experiment acted in a comparable manner in the selected stands. It is likely that these large-scale factors act in a more or less unequal manner in each selected stand, which could have some influence on oak regeneration density. To address this issue, mixed effects models were used in the statistical analysis to separate the variability of the dependent variable attributable to the factors studied (fixed effect, e.g. local stand basal area) from the variability caused by factors not directly included in the model that are specific to the particular forest stand (random effect).
Starting from a random point, we established sample plots in a rectangular grid fitted to the area of each stand. The total number of sample plots in all stands was 240 (on average 18.5 per stand). The measurements were carried out on three concentric circular plots of different radii. Seedlings of oak and other tree species (h ≤ 0.5 m) were counted only in the smallest area of 10 m 2 (radius 1.78 m). On the second circle (area of 100 m 2 , radius 5.64 m) saplings (h > 0.5 and diameter at breast height (DBH) ≤ 2 cm) and recruits (DBH: 2-7 cm) were counted. Trees (specimens with DBH > 7 cm) were measured in the largest circle (area of 250 m 2 , radius 8.92 m). We measured the DBH of all trees, the height of trees, and the radii of the tree crown in four directions using tape then, the crown projection area was calculated using formula for ellipse area. In the case of regeneration, we measured the height of all seedlings and saplings and the DBH of saplings taller than 1.3 m. The quality of the recruit stem was categorized into the following classes: 1straight, 2curved (one or more curves), and 3deformed (bushy shape). In further analysis we divided recruits into two quality classes: first quality (straight) and combined low quality (curved and deformed). The cover (%) of litter, moss, and herbs was described for each circular plot (in the area of 250 m 2 ).
Soil enzyme activity data were collected in 10 sample plots (the same plots established for stand measurements) in each stand studied. After removing the litter layer (Ol layer), soil samples were taken and collected from the top 15 cm of soil (this is the depth where the main microbial activity occurs) using a 200 cm 3 probe. The soil samples were divided into two soil horizons: organic (Oh holorganic layer) and humus (Ah organomineral horizon). The soils were then sieved (< 2 mm), removing all debris, stones, roots and plant remnants. Samples taken from one stand were pooled and further processed to obtain an estimate representing the entire stand.
Samples taken from one stand were pulled together and processed further to receive one estimate representing the whole stand. Four enzymatic activities were analysed (dehydrogenase, urease, phosphatase and asparginase). Dehydrogenase activity was determined by the reduction of 2,3,5-triphenyltetrazolium chloride (TTC) to triphenyl formazan (TPF) using Lenhard's method according to the Casida procedure (Alef and Nannipieri 1995). Briefly, 1 g of soil was incubated with 1 mL of 3% TTC for 24 h at 37°C. TPF was extracted with ethyl alcohol and measured spectrophotometrically. Urease activity was determined according to Tabatabai and Bremner (1972) using a water-urea solution as a substrate. This activity was determined by the NH 4 + released after a 48 h incubation at 37°C. The concentration of NH 4 + was measured at 410 nm by the colorimetric method (Alef and Nannipieri 1995). The asparginase activity was determined with the use of the colorimetric technique and expressed as NH 4 + mg in 10 g of soil per 48 h. The acid phosphatase activity was determined by the calorimetric method and expressed in mg of p-nitro phenol (pNP) per 10 g of soil (Tabatabai and Bremner 1969;Olszowska 2018).

Statistical analyses
Generalized linear mixed models with a log-link function were developed for five separate response variables (oak seedling counts, oak sapling counts, oak recruit counts first grade oak recruits counts and summarised medium and lower grade recruits counts) to assess factors influencing oak numbers belonging to particular regeneration stages observed on the sampling plot. Although sample plots were placed on systematic grids, plots placed in one forest stand could be more related to each other (dependent) than to other plots. Additionally, some of the independent variables (e.g., describing soil enzymatic activity or the age of Scots pine trees on sample plots) were estimated at the stand level. To disentangle potential interactions between many factors and address autocorrelation, generalized linear mixed effect models (GLMMs) were used (Zuur et al. 2009;Bolker et al. 2009). The influence of the common location of a group of sampling plots from one forest stand was modelled as a random effect, and other potential factors were modelled as fixed effects. In the mixed models, the intercept for the random effect was allowed to vary between forest stands, but the slope was fixed.
Count data are discrete and nonnegative; thus, the Poisson distribution is often used for modelling such data. An important limitation of the Poisson distribution is the assumption that the mean and variance of the observed data are equal. The situation that the observed variance is greater than the mean is quite common in ecological data sets and is referred as over-dispersion (Zuur et al. 2009). To account for the possibility of such a situation in oak regeneration modelling, a negative binomial distribution of the Poison distribution (GLMM Poisson) was also used. The negative binomial distribution can be parameterized differently regarding the dependence of the variance on the mean (Crotteau et al. 2014). In the present study, two variants of dependence were implemented: one when the variance increases linearly with the mean (GLMM nbinom1) and the second when the variance increases quadratically with the mean (GLMM nbinom2).
Some histograms from Fig. 1 suggest that the proportion of zeroes is so large that it could not be possible to readily fit these data to standard (Poisson or binomial) distributions. The non-consideration of this excess of zeroes might reduce the ability to detect relevant relationships and make inaccurate inferences. Three additional types of zero-inflated mixture models were built: one based on a Poisson distribution (GLMM ZI Poisson) and two based on a negative binomial distribution with linear (GLMM ZI nbinom1) and quadratic (GLMM ZI nbinom2) relations between mean and variance. Finally for each regeneration stage, six additional mixed models for each oak regeneration stage were built (GLMM Poisson, GLMM nbinom1, GLMM nbinom2, GLMM ZI Poisson, GLMM ZI nbinom1, and GLMM ZI nbinom2).
All the possible explanatory variables collected during field work are shown in Additional file 1: Table A1. In the first stage of modelling simple generalized linear Poisson models containing only one explanatory variable were build. About twenty variables which gave the best models were selected to the next modelling stage in which all of them were included in a model belonging to one of the previously mentioned model types (e.g. GLMM nbinom1). In the case of collinearity, the variable with less explanatory power were removed from initial models based on their inflation factors (VIFs) (Zuur et al. 2009). In the refined set of explanatory variables no variable had VIFs value above 3. During building of statistical models interaction terms were added to assess potential influence of stand characteristics on the reactions of young oaks on the change of soil enzyme activity level. For building models concerning recruit counts, variables concerning vegetation soil coverage were not included because it is rather unlikely that the number of trees with a DBH greater than 2 cm is governed by low vegetation.
At the end of the model fitting process for each oak regeneration stage, six competing models were built that address the problems of autocorrelation, overdispersion and zero inflation. The final model describing the influence of important variables on the particular oak regeneration stage was chosen based on the lowest AIC value from 6 models (Zuur et al. 2009(Zuur et al. , 2012Crotteau et al. 2014;Peters and Visscher 2019). After the identification of the best model, its residuals were inspected to confirm that errors were homogeneous.
All analyses were performed with R language (R version 3.5.0, R Core Team 2018), and models were fitted using a template model builder (TMB) via maximum likelihood estimation using the R package "glmmTMB" (Brooks et al. 2017). For the final models containing random effects, it was possible to characterize their explanatory power by the calculation of the marginal R-squared value considering only the variance of the fixed effects and conditional R-squared, taking both the fixed and the random effects into account according to the formulas proposed by (Nakagawa et al. 2017).The R 2 values for the selected models were calculated with the function from the "performance" package (Lüdecke et al. 2019).

Results
The general information of specimen counts belonging to different species and developmental stages are given in the Additional file 1 in Tables A3-A5. The frequencies of different classes of oak regeneration on sampling plots and information on the overall density of oak spontaneous regeneration are presented in Fig. 1. We found that the density of oak regeneration decreased from the seedling class to larger classes. The total density of oak seedlings was almost two times higher than the sapling density (2043 vs. 801 individuals per ha). The number of recruits was lower than 200 individuals per ha. We compared the density of first quality (straight stems) and lower quality (curved stems and bushy shape) recruits and found that the density of lower quality recruits was greater than that of better quality recruits.
In this section, we refer only to the final models. Competing model details and comparisons can be found in the Additional file 1 (Tables A6-A10). The assumption that spatial autocorrelation connected with plot cooccurrence in the same stand was important was supported by lower AIC values of models containing stand random effects for almost all final models but not for the best model describing seedling counts. The second best competing model for seedlings was the GLMM with a comparable log likelihood value, but its AIC value was higher by 2 due to a larger degrees of freedom value, which means that including random effects in this case did not improve the model. In no final models interaction term between soil enzyme activity and other independent variables was proven to be statistically important.
The usefulness of the Poisson distribution for modelling oak regeneration was very low. In none of the analysed regeneration stages were competing models based on the Poisson distribution found among the three highest rated models (Additional file 1: Tables A6-A10). In the majority of cases, Poisson models occupied the end of the ranking table. The over-dispersion value for the best models ranged from 1.29 (in the recruit model) to 3.29 (in the sapling model). Three of the five best models were built on the assumption that the dependence between the mean and variance of the specimen count had a nonlinear form (nbinom2), but the examination of competing model ranking tables (Additional file 1: Tables A6-A10) suggested that there was no unequivocal dominance of this distribution form among the highest ranked models.
Only for seedling model the assumption concerning zero inflation was proven by the AIC comparison of competing models. The best mixed model build for seedlings counts namely GLMM ZI nbinom2 turned out to have variance of random component equal zero. This problem was solved in the way described by Pasch et al. (2013) by the reduction of random component of the model. The parameters of the final model were recalculated back from logit form and provided in Table 1. The only explanatory variable proven to be useful for predicting zero excess was bilberry (Vaccinium myrtillus L.) ground cover. The ground cover of raspberry (Rubus idaeus L.) growing on the sample plot and the birch crown projection area had a positive influence on the observed number of oak seedlings. A negative effect was cast by the ground cover with litter and ferns (Pteridium aquilinum (L.) Kuhn).
From analysed soil enzymes only activity of dehydrogenase correlated positively and in statistically important manner with numbers of seedlings on sampling plots (Fig. 2b). If the dehydrogenase activity increased to 0.4, the predicted average count of oak seedlings was greater than 3. The increase in the raspberry cover to 25% and birch crown projection area to 105 m 2 (the highest observed levels of those factors) strongly positively influenced oak seedling counts (Fig. 2a, c). The negative impact of fern and litter was presented on Fig. 2e, d. Plots without litter and fern cover had markedly higher seedling counts.
Four predictors had a negative influence on the number of oak saplings (Table 2). There were: ground cover by litter, basal area of pedunculate oak, European hornbeam and Norway spruce. A positive relationship was observed for phosphatase activity in the organic soil horizon. As shown in Figs. 3a-d the negative influencing factors have a very strong diminishing influence on the saplings counts. The highest oak sapling density was predicted if the basal area of oak, spruce and hornbeam was zero. The higher crown projection area of these species decreased the oak sapling count. Compering the effect of these species, the strongest impact was observed for hornbeam. We found that if the litter cover was 100% the count of oak saplings was 0.6. As shown in Fig. 3e if the phosphatase activity is greater than 4.5 the expected count of oak was 7.
Three statistical models describing counts of different types of oak recruit were built. The first (Table 3) described general counts of recruits. This model contained only two statistically important variables (the activity of phosphatase in the humus soil horizon and total crown projection area of all trees growing on sample plots) that negatively reacted with the predicted recruit counts. The predicted count of all recruits were presented in Fig. 4a and b. If the total crown area was 400 m 2 the count of oaks was close to zero. As shown in Fig. 4b, only when the phosphatase activity in the soil was close to the lowest observed levels (0.23 mg per 10 g of soil) all recruit counts was 2.5. A similar pattern could be observed when only the first quality (in terms of their potential silvicultural use) recruits were analysed (Table 4, Fig. 4c and d). The same factors influenced recruit counts in a similar manner, but the average count of this class of    recruits was two times lower than that of all recruits. In the case of the lower quality recruit counts (Table 5, Figs. 4e-g), a statistically important negative impact of three predictors was detected: the activity of phosphatase in the humus soil horizon, spruce basal area and oak volume in the sample plot.

Discussion
Today's challenge in forestry goes further than to achieve timber production with social requirements and to ensure maintenance of forests as part of our heritage (Schütz 1999; O'Hara 2016). Expected climate change will influence vast forest areas simultaneously and in a relatively short time in comparison to the longevity of the classical forest production cycle. The gradual adjustment of species composition during the process of utilisation of mature forest stands and establishment of new better adapted generation could be too slow to react on fast occurring climate changes (Schelhaas et al. 2015). Near-natural silviculture approaches emphasizing the utilization of natural spontaneous processes may have the potential to develop complex and sustainable forests that are adapted to our changing world (Brang et al. 2014). The idea that biological processes rather than silvicultural efforts might be relied upon are especially promising given the size of the challenges ahead. Biological automation might accommodate reduced resource inputs/effort into forest management operations (Pretzsch and Zenner 2017). Our investigation referred to the idea of creating mixed stands that are well adapted to changing environments. Moreover, most studies have focused on the most popular 2-species combinations (e.g., spruce-beech, oak-beech), while other important combinations, such as pine-oak, have received scant attention ). Because climate change may negatively impact growing conditions for Scots pine monocultures situated on dry, sandy soils in central Europe (Slodicak et al. 2011), spontaneous oak regeneration can play an important role in their transformation to future stands that are more stable and species reach.

Density of oak regeneration
We observed oak regeneration of all development phases (from seedlings to recruits) in Scots pine-dominated stands (Additional file 1: Tables A3-A5). Half of all oaks were included in the smallest height group (up to 0.5 m height). We would like to stress that in the selection of important variables in the model building process, we did not show that the age of Scots pine stands is an important variable.
The average density of oak recruits was not high, but the importance of their presence was potentially great, especially if they were present in the relatively young Scots pine-dominated stands (Additional file 1: Table  A5; Fig. 1). Although their direct economic meaning is rather constrained (they could be utilized at most for fuel wood during final cuttings), their impact on the biological stabilization of Scots pine stands is well documented (Steckel et al. 2019). The conversion of pure Scots pine into mixed stands with understory oaks soundly increased the number and species composition of parasitoid wasps, which could mitigate the outbreaks of folivore insects (Jäkel and Roth 2004). This goal could be achieved in artificial way by underplanting oaks in 40-50 year old Scots pine stands which is costly, or in natural way by jays activity which we demonstrated in our results.
Some modelling results suggested that even recruits that were not abundant could potentially have economic importance because they are not randomly but rather contagiously distributed in the stand. In the beginning stage of final model selection, initial statistical models assuming different types of individual distributions (e.g., Poisson, negative binomial and their zero inflated counterparts) were compared for each regeneration stage  Tables A6-A10). As in many other investigations (Fyllas et al. 2008;Li et al. 2011;Zhang et al. 2012;Crotteau et al. 2014;Vickers et al. 2017), models that did not assume Poisson distributions of young regenerating individual counts proved to be better. Although the better performance of models based on a negative binomial count distribution is not a general pattern in analysing spontaneous regeneration (Vickers and Palmer 2000;Fyllas et al. 2008;Peters and Visscher 2019), this was an important finding from the management (silviculture) point of view. The Poisson distribution of counts on randomly chosen sampling quadrats could be seen as evidence that the spatial point pattern formed by specimens is completely random (Diggle 2003;Wiegand and Moloney 2014). The eventual superiority of Poisson-based models describing oak regeneration counts on sampling plots would suggest that factors influencing oak regeneration spatial placement were acting in a random manner on the area of investigated stands. The superiority of models based on different forms of negative binomial distributions suggested that regenerating oaks were not placed completely randomly. Based only on this information, it was not possible to note the name of a particular spatial point process "responsible" for creating oak spatial distribution (Diggle and Milne 1983;Coly et al. 2016), but it was safe to assume that the spatial distributions of oak were heterogeneous (Velázquez et al. 2015) or even grouped (Krebs 1999) in the space of investigated stands. Grouped distribution of different stages of birddispersed oak regeneration was also found in other studies (Mosandl and Kleinert 1998;Frost and Rydin 2000). From a practical point of view, this result means that locally (in some fragments of a forest stand), the density of oak saplings or recruits could be large enough to be useful for silvicultural goals, e.g., in some forest districts in Poland, clumps of well-shaped oak recruits are used to create the oak admixture in the next generation of Scots pine stands (Gniot 2007;Skrzyszewski and Pach 2015). The possibility of successfully incorporating understory oaks as a good quality admixture into the next generation is limited by the amount of time young oaks spend in the understory. The longer the time, the greater the chance that their stems will become crooked. When planning the conversion of Scots pine to oak, full overstory light should be provided as early as possible, but no later than 20 years after the regeneration is established (Skrzyszewski and Pach 2015).

Soil factors influencing oak spontaneous regeneration
We found that enzymatic activity is one of the most important factors correlating with oak spontaneous regeneration in Scots pine stands. Our statistical models suggest that the correlation with soil enzyme activity depends on the stage of oak regeneration. The number of oak seedlings was positively related to the soil dehydrogenase activity. This is often interpreted as an indicator of increased microbial activity in the soil, particularly mycorrhizal fungal activity (Buée et al. 2005). Yguel et al. (2014) found that oaks surrounded by phylogenetically distant neighbours had increased abundance and enzymatic activity of ectomycorrhizal fungi in the litter. This suggests that reduced nutrient availability in a phylogenetically distant litter was partially compensated  Observations 240 Marginal R 2 / Conditional R 2 0.500 / 0.723 for by increased litter decomposition by ectomycorrhizal fungal activity. The research conducted by Showalter et al. (2010) directly shows that dehydrogenase activity is positively correlated with tree growth, which, according to the cited author, indicates that a well-established mycorrhiza is increasing nutrient availability for host tree.
The density of oak saplings and recruits was related to the phosphatase activity in soils. The increase in phosphatase activity in organic soil horizon corresponded to increased oak sapling density but the density of oak recruits was negatively correlated with phosphatase activity in the humus soil horizon.
Soil phosphatase plays critical roles in phosphorus cycles and the metabolic state of soil microorganisms (Watts et al. 2010) and its activity is positively correlated with soil-extractable phosphorus and with high production capacity, stand biomass and/or plant cover (Carreira et al. 2000). Phosphorus availability is essential for plant growth and may be a limiting factor in some forest ecosystems (Attiwill and Adams 1993). This constraint could be especially important for young oaks (Collet et al. 1997), with small and relatively shallow root system growing in relative poor site condition in our experiment. The higher activity of the enzyme increases the nutrient uptakes from organic soil horizon by oak regeneration and could promote their survival and increase the density of oak saplings.
We hypothesise that young oaks depend positively on soil enzyme activity, but for the older ones the causeeffect relationship is reversed, so that soil enzyme activity depends negatively on the presence of recruits. The size of the recruits is much larger than that of the other oak regeneration classes studied. These relatively large organisms could have a more significant effect on soil microbial activity than smaller ones. There are at least two ways in which oak recruits may reduce phosphatase activity in Scots pine-dominated stands: 1) Oak trees took up phosphorus mainly from 15 cm soil depth, where the greater amount of roots and external mycorrhizal mycelia were found (Göransson et al. 2006). The phosphatase activity is correlated with the availability of phosphorus in the soil. The increase in phosphorus in soil typically leads to a decrease in the activity of this enzyme (Olander and Vitousek 2000). The additional amount of available phosphorus in soil coming from the decomposition of litter composed not only with pine needles but also oak leaves could decreases the activity of phosphatase. Soils under Quercus typically showed low enzyme activity (Šnajdr et al. 2013). Although the amount of phosphorus is similar in pine and oak litter and soil (Šnajdr et al. 2013), the rate of their litter decomposition is different. Oak litter is rapidly decomposed compared to other litter (Šnajdr et al. 2013). Moreover, the mineralization rate is higher in soils under broadleaved trees than under Norway spruce and Scots pine Smolander and Kitunen (2002). For oak younger regeneration growing in soil covered by almost pure Scots pine litter the higher activity of phosphatase is needed to improve the supply of phosphorus. Oak recruits add substantial quantity of leaves to the litter so availability of phosphorus could be improved by the higher ratio of litter decomposition and the increased activity of phosphatase is not needed.
2) The second mechanism in which the presence of oak recruits could negatively influence the activity of phosphatase could be connected with their influence on soil moisture. Oaks have different strategy than Scots pine in terms of soil water usage. In the case of drought Scots pines strongly reduce their transpiration but oaks in Central Europe tend to keep high rate of transpiration as long as possible (Toïgo et al. 2015). This difference in water usage strategy was observed in Europe-wide experiments (Steckel et al. 2020) and could lead to faster water depletion under oaks than pines. Augusto et al. (2003) observed that vascular plants growing under oaks in mixed stands have lower moisture requirements than under Scots pine. The presence of oak recruits in investigated stands could locally diminish soil moisture and diminish activity of phosphatase, which depends strongly on this soil property (Baldrian 2014).

Impact of stand features on the density of oak regeneration Oak seedlings establishment
We found that the cover of bilberry negatively influenced the establishment of oak seedlings. However, Drössler et al. (2017) observed more oaks in blueberry patches and suggested that the Eurasian jay (Garrulus glandarius) prefers to hide acorns under dwarf shrub vegetation. A negative impact on oak seedlings was observed from fern cover. Jensen et al. (2011) also suggested a negative effect of dense herbaceous ground vegetation on oak regeneration. Bilberry and fern, especially common bracken, created dense ground cover that was not preferred by oak seedlings in pine stands. Humphrey and Swaine (1997) showed that competition from bracken (Pteridium aquilinum (L.) Kuhn) restricted the growth of oak seedlings. Competition for nutrients and moisture may also be important, especially in nutrient poor or drier areas (Löf 2000;Brudvig and Asbjornsen 2007). However, we also found that litter cover created inappropriate conditions for oak regeneration. In previous studies (Kurek and Dobrowolska 2016;Kurek et al. 2018), it was observed that jays deposited the acorns in small patches of the litter. Moreover, litter affected soil humidity and the amplitude of diurnal temperature fluctuations. The forest floor can act as a physical barrier and can possibly release toxic metabolites (Sayer 2006). This process depends on the amount of litter present and the environmental conditions (Donath and Eckstein 2008). The positive impact of raspberry cover on oak seedling density was observed in our study. Raspberry does not create as dense of ground cover as bilberry, and it protects oaks against damage caused by biotic (herbivory) and abiotic (drought, insolation, lack of humidity) factors (Donoso and Nyland 2006). Cuttings that remove more than 40% of the forest canopy create environmental conditions that promote the establishment of raspberry. In such places, oak finds good conditions for regeneration (Donoso and Nyland 2006).
Overstory species composition affected oak seedling density. We recognized the positive impact of birch (birch crown projection area) on oak density. Light transmission was found to be higher in dense birchdominated stands than in dense pine-dominated stands because of the higher total foliage area and the higher location of foliage in the pine canopy (Lintunen et al. 2013). Because pedunculate oak is a lightdemanding species (Savill 2019), it requires at least 20% full sunlight to avoid severe growth depression (Ligot et al. 2013). Light is not a requirement for germination (Ligot et al. 2013), as seedlings largely rely on energy from the acorn during the first season. Paluch and Bartkowicz (2004) also found that oaks occurred more frequently in the vicinity of birches. It is possible that the neighbourhood of birch trees could facilitate the establishment of oak by reducing the competition of vegetation.

Oak saplings
We recognized that hornbeam, Norway spruce and pedunculate oak basal area negatively influenced the number of oak saplings. All of these species have dense crowns that transmit less light. We think that light conditions were the key factor, as the light requirement of oak increases with increasing tree age and size (von Lüpke and Hauskeller-Bullerjahn 1999;Vizoso-Arribe et al. 2014). Annighöfer et al. (2015) also showed that the occurrence of oak saplings was related to light conditions and that abundance increased with increasing light availability. Moreover, the negative impact of litter cover on sapling density suggested more demands of advanced oak regeneration for light. Our results were opposite to those of Lithuanian investigators, who observed that the abundance of oak undergrowth was largest where spruces and beeches were predominant in the overstory (Jurkšienė and Baliuckas 2018). However, their results confirmed our results regarding the negative influence of hornbeam on oak regeneration (Jurkšienė and Baliuckas 2018). Although pedunculate oak can grow as a pioneer tree (Götmark and Kiffer 2014), it survives as a seedling/sapling in relatively dark understories (our results).

Oak recruits
The count of oak recruits showed a clear relation to crown projection of all trees. The density of oak recruits increased with decreasing crown projection of all trees. Crown size is an indicator of space occupancy because it is correlated with the photosynthetic capacity (Kamler et al. 2016). Similar results were achieved by Annighöfer et al. (2015), who found that sapling quantity decreased with increasing basal area of other species. The situation when the number of recruits is low and the crown projection of old tree is high suggests water limitations at small scales. On the Fig. 4, it could be seen that at the same level of local crown projection, the density of high quality recruits was lower than the total density of recruits. Recruits is the category of young oaks that have the greatest age and longer growth history under canopy pressure. Observation from other studies (Skrzyszewski and Pach 2015) indicate that prolonged period of growth under the canopy (more than 20 years) reduces the quality of young oaks. This explanation is suitable also for the first-quality recruits. In the case of low-quality recruits, more factors had a negative impact on their density. These factors included spruce basal area and oak volume. Norway spruce and oak utilise more water than Scots pine and may locally diminish water reserves in soil and also transmit less light through their canopy. With the lack of light, oaks create a shrubby crown. Tall and slender oaks reflect a priority for shoot growth, which is a common strategy employed by plants in response to shading (Jensen et al. 2011).

New dimension of oak regeneration niche
In our study, we explored the oak regeneration niche, i.e., the set of environmental requirements potentially important for germination and establishment of its regeneration. Much research is devoted to exploring variables that constrain oak regeneration, and indeed they explore various dimensions of the oak regeneration niche, but rarely do authors directly state that they are studying this phenomenon (Collins and Good 1987). The concept of the regeneration niche (Grubb 1977) is based on the idea of the ecological niche, which was defined by Hutchinson (1957) as a region in a multidimensional space of environmental factors that influence the well-being of a species. The results of our study reveal new, potentially important dimensions of the oak regeneration niche. The relatively high coefficients of determination of models describing the number of young oaks in the stands studied suggest that soil conditions, represented by soil enzyme activity, play an important role in the establishment and growth of oak spontaneous regeneration under the canopy of Scots pine stands. This factor, relatively rarely studied in the context of oak regeneration, may be important in explaining the failure of oak regeneration in some central European forests. It is well known that the importance of a particular dimension of ecological niche may change in different parts of a species' distribution range, e.g. in a colder region the availability of direct sunlight might be more important than in a warmer one (Peterson et al. 2011). It is likely that the importance of soil enzyme activity for oak regeneration establishment could change with changes in other environmental variables, especially beyond the boundaries examined in our study. It is difficult to give simple instructions on how to regulate the level of enzyme activity on an economic scale. However, we believe that the results of our study can be used to some extent to diagnose the proper site conditions for oak regeneration.

Conclusions
Successful regeneration of pedunculate oak under Scots pine-dominated stands of different ages is possible, as shown by the presence of all stages of oak regeneration from seedling to recruit. We found that oak regeneration density depended on a combination of several variables, but the activity of two soil enzymes played a major role in oak establishment and advancement. Soil enzyme activity can be considered not only a predictor of site conditions, but also a predictor of establishment and advancement of oak regeneration. The results of our study reveal new, potentially important dimensions of the oak regeneration niche.
The spatial distribution of oak saplings and recruits was heterogeneous or even grouped in some fragments of a studied forest stand; therefore, it can be used in future conversion of pine-dominated stands into mixed stands. Even if oak cannot be considered an important tree species in the upper stand layer now, it will play an important role in the future forest ecosystem because existing groups of good quality oaks (especially saplings) could be used as admixture when creating the next generation of forest stands. It could be particularly useful if possible climate change forces us to convert large areas of Scots pine monocultures in Central Europe into mixed forest stands. Even if it will not be necessary, the spontaneous spread of oaks in Scots pines monocultures could increase the biological stability and resilience of these forests (resistance to outbreaks of folivorous insects).
Abbreviations DBH: Diameter at breast height; TTC: 2,3,5-triphenyltetrazolium chloride; TPF: Triphenyl formazan; pNP: p-nitro phenol; GLMMs: Generalized linear mixed effect models; GLMM Poisson: Negative binomial distribution of the Poison distribution; GLMM ZI Poisson: Zero-inflated mixture model based on a Poisson distribution; GLMM nbinom1: Mixture model based on a negative binomial distribution with linear relations between mean and variance; GLMM ZI nbinom1: Zero-inflated mixture model based on a negative binomial distribution with linear relations between mean and variance; GLMM nbinom2: Mixture model based on a negative binomial distribution with quadratic relations between mean and variance; GLMM ZI nbinom2: Zero-inflated mixture model based on a negative binomial distribution with quadratic relations between mean and variance
Additional file 1: Table A1. Expalanatory variables collected during field work. Table A2. Species composition of the investigated stands. Table A3. Density of seedlings (H ≤ 0.5 m) in the investigated pine dominated stands. Table A4. Density of saplings (H > 0.5 m and DBH ≤ 2 cm) in the investigated pine dominated stands. Table A5. Density of recruits (DBH < 7 cm) in the investigated pine dominated stands. Table  A6. Comparison of competing models for seedlings on Akaike Information Criterion (AIC). GLM = generalized linear model; GLMM = generalized linear mixed model; ZI = zero inflated; Family = family of error distribution; negbin1 = negative binomial (variance that increases linearly with the mean); negbin2 = negative binomial (variance that increases quadratically with the mean); Df = degrees of freedom; The optimal model is placed on the top of table. dLogLik and dAIC are the difference between subsequent models and the best one in term of AIC and log likelihood (logLik). Table A7. Comparison of competing models for saplings on Akaike Information Criterion (AIC). GLM = generalized linear model; GLMM = generalized linear mixed model; ZI = zero inflated; Family = family of error distribution; negbin1 = negative binomial (variance that increases linearly with the mean); negbin2 = negative binomial (variance that increases quadratically with the mean); Df = degrees of freedom; The optimal model is placed on the top of table. dLogLik and dAIC are the difference between subsequent models and the best one in term of AIC and log likelihood (logLik). Table A8. Comparison of competing models for recruits on Akaike Information Criterion (AIC). GLM = generalized linear model; GLMM = generalized linear mixed model; ZI = zero inflated; Family = family of error distribution; negbin1 = negative binomial (variance that increases linearly with the mean); negbin2 = negative binomial (variance that increases quadratically with the mean); Df = degrees of freedom; NA = not applicable to computational/convergence issues. The optimal model is placed on the top of table. dLogLik and dAIC are the difference between subsequent models and the best one in term of AIC and log likelihood (logLik). Table A9. Comparison of competing models for first quality recruits on Akaike Information Criterion (AIC). GLM = generalized linear model; GLMM = generalized linear mixed model; ZI = zero inflated; Family = family of error distribution; negbin1 = negative binomial (variance that increases linearly with the mean); negbin2 = negative binomial (variance that increases quadratically with the mean); Df = degrees of freedom; NA = not applicable to computational/convergence issues. The optimal model is placed on the top of table. dLogLik and dAIC are the difference between subsequent models and the best one in term of AIC and log likelihood (logLik). Table A10. Comparison of competing models for lower quality recruits on Akaike Information Criterion (AIC). GLM = generalized linear model; GLMM = generalized linear mixed model; ZI = zero inflated; Family = family of error distribution; negbin1 = negative binomial (variance that increases linearly with the mean); negbin2 = negative binomial (variance that increases quadratically with the mean); Df = degrees of freedom; NA = not applicable to computational/ convergence issues. The optimal model is placed on the top of table. dLogLik and dAIC are the difference between subsequent models and the best one in term of AIC and log likelihood (logLik).