Seedling damage caused by wood harvesting and soil scarification in rotation and continuous cover forestry in Scots-pine-dominated boreal forests

ABSTRACT Continuous cover forestry (CCF) is gaining popularity as an alternative to rotation forest management (RFM), especially in forests with multiple uses. The success of CCF depends on the amount and quality of remaining seedlings, but the effect of harvesting on their status is not well known. In the present study, the proportion and number of undamaged seedlings after harvesting (and soil scarification) in Scots-pine-dominated stands in Northern Finland were modelled by applying logistic binomial mixed-effects and negative binomial count data models. The treatments included both CCF (gap cutting and selection cutting) and RFM (clearcutting and seed tree cutting) harvesting methods. The clearcut areas and most seed tree areas were treated with disc trenching. In the clearcutting and seed tree cutting areas, 25–34% of the inventoried seedlings were damaged. In the gap cutting and selection cuttings, c. 7% and 12% of the seedlings were damaged respectively. Harvesting machinery traffic and soil scarification were the most important causal agents of seedling damage. An increase in slash coverage also significantly increased seedling damage. The average number of undamaged seedlings in the selection cutting areas was ca. 1500 pcs ha−1, indicating at least satisfactory regeneration potential in the experimental stands.


Introduction
Continuous cover forestry (CCF) as an alternative to rotation forest management (RFM) based on final fellings is driven by environmental concerns, national and EU legislation, forest certification programmes, and various EU strategies (Laudon and Maher Hasselquist 2023).CCF includes silvicultural systems involving continuous and uninterrupted maintenance of forest cover without clearcutting (Pommerening and Murphy 2004).CCF has been demonstrated to maintain biodiversity (Peura et al. 2018) and carbon stocks (Assmuth and Tahvonen 2018), though there are also studies that have not observed such effects (Koivula 2012;Lundmark et al. 2016;Rikkonen et al. 2023).Avoiding clearcuts and soil scarification improves the quality of runoff waters (Laudon and Maher Hasselquist 2023).In peatland forests, CCF could be an alternative approach for managing groundwater levels by maintaining the evapotranspiration capacity of the trees without ditch network maintenance.Moreover, greenhouse gas emissions from peat soils could be counteracted by regulating groundwater levels with CCF (Nieminen et al. 2018).Only 2.7% of the forests in Finland in 2020 were managed with CCF or were undergoing transformation to CCF.The proportion of the area harvested with CCF methods was highest (5.8%) in Lapland (Finnish Forest Centre 2020).
In addition to the ecosystem services listed above, CCF could also facilitate the integration of wood production with multiple uses of forests (Miina et al. 2020).For example, tourism, reindeer herding, and forestry are important livelihoods in northern Fennoscandia.Parkatti and Tahvonen (2021) have shown that including the negative effects of forestry on reindeer husbandry in economic calculations favours continuous cover forestry in northern Scots pine stands.The soil disturbance caused by harvesting and soil scarification machinery reduces the coverage of ground lichens, which are an essential element of reindeer's winter nutrition (Harris 1992;Roturier and Bergsten 2006;Kumpula et al. 2008;Kivinen et al. 2011).Moreover, the slash (harvesting residue) accumulating on the ground has been observed to suppress lichens (Bråkenhielm and Liu 1998;Kumpula et al. 2008;Akujärvi et al. 2014).Consequently, cautious harvesting of mature forests typical of CCF could better maintain goodquality winter pasturelands than RFM does (Hirvelä et al. 2022).Continuous cover logging methods have also been recommended on sites close to settlements and recreational areas, as they maintain the aesthetic and recreational values important for nature-based tourism (Koivula et al. 2020).
Stand structure in CCF can be uneven-aged, or the trees can be grown in two or more storeys.In Finland, selection and gap cutting, shelter tree cutting, thinning of hold-overs, opening-ups, and shade tree cutting are the main cutting methods associated with CCF.The partial cuttings increasingly applied in North America are analogous to CCF, as their aim is to maintain age structure and tree diameter distribution in accordance with the principles of ecosystem-based management and to increase tree growth by reducing competition (Caspersen 2006;Pamerleau-Couture et al. 2015).Cuttings promoting the transformation from RFM to CCF are also included in the CCF practices (Valkonen 2022).Selection cutting is defined as the annual or periodic removal of trees (particularly mature ones), individually or in small groups, from an uneven-aged forest to achieve the balance among the diameter classes required for a sustained yield and to realise the yield and establish a new crop of an irregular constitution (Nieuwenhuis 2000).Gap cutting is used to promote the emergence and growth of seedlings.In Finland, the maximum area of an individual gap is limited to 0.3 ha (Forest Act 2014), but gaps can be enlarged at subsequent harvest times after successful regeneration in gaps.
The economic competitiveness of CCF compared to RFM is enhanced with increasing discount rates and management and with decreasing site productivity and timber prices (Tahvonen 2009).In extreme areas such as Lapland, artificial forest regeneration based on planting or seeding with subsequent tending can impair the cost-effectiveness of wood production.In CCF wood harvesting greatly affects the future development of uneven-aged stands, as the stands are expected to recover by the next harvesting phase through the growth of remaining stand and stand regeneration facilitated by the removal of trees (Surakka and Sirén 2007).In RFM in Finland stands are typically regenerated by clearcutting with subsequent soil scarification and cultivation or by seed tree cutting.Clearcut areas can be planted or seeded; seeding is usually integrated with soil scarification.In seed tree cutting natural seeding by the seed trees is expected.However, in northern areas good seed crops of Scots pine (Pinus sylvestris L.) and Norway spruce (Picea abies L. (Karst)) are rare due to the cool climate (Juntunen and Neuvonen 2006;Hilli et al. 2008;Pukkala et al. 2010).Birch (Betula pendula Roth and Betula pubescens Ehrh.) are prolific seed producers, but there is a large annual variation in the quality and quantity of the seed crop (Hynynen et al. 2010).Avoiding harvesting damage to the undergrowth is therefore very important, particularly in northern regions, where climate conditions limit natural forest regeneration.The risk of biotic and abiotic damage to forests is also increasing due to climate change (Venäläinen et al. 2020).
In Mason et al.'s (2022) survey, the use of mechanised harvesting systems in CCF was ranked as one of the main knowledge gaps limiting the implementation of CCF in countries where CCF had only recently been accepted, and the forestry was based on even-aged stand management through RFM (e.g.Finland).The cut-to-length method (CTL) is the dominant harvesting method in Northern Europe.Lundbäck et al. (2021) estimated that mechanised harvesting accounts for up to 95% of the total harvest in Finland, Sweden, and Norway.Trees are usually cut with a single-grip harvester that fells and crosscuts the stems into timber assortments on site.In most cases, the timber is transported to the roadside using forwarders.Part of the productive area is lost because of the strip roads required for machine access routes and timber transport.According to the silvicultural guidelines for RFM in Finland, the minimum allowable strip road spacing in intermediate cuttings is 19 m, and the allowable width for strip roads 4.6-5.1 m (Leivo et al. 2022).In the final fellings, strip road spacing varies from 10 to 15 m (Nurminen et al. 2006).
Harvesting techniques and machine constructions are developed for even-aged forestry, in which saving the undergrowth is not usually required (Surakka and Sirén 2007).Understorey trees hinder visibility when selecting the trees to be cut.They can also cause technical problems for the harvester (e.g.breakages of saw chains and hydraulic hoses) and prevent the operator correctly positioning the harvester head (Kärhä and Bergström 2020).Due to the decline in the efficiency of operations, undergrowth trees obstructing harvesting are typically cleared before cutting in RFM.
In general, an increase in removal volume increases damage to the remaining trees (Nuutinen and Muhonen 2022).Felling causes more than 60% of the damage to the remaining trees, and the proportion of delimbing and crosscutting exceeds 20% (Sirén 1998).The felling of large trees above the undergrowth in selection cutting causes stem and top breakages.Moving the felled trees back and forth in the harvester head for delimbing and crosscutting damages the undergrowth.The processing of large birches with wide tops especially presents a risk to the remaining stand.The piling of logging residues and roundwood on the undergrowth further damages the undergrowth trees (Surakka and Sirén 2007;Surakka et al. 2011).In North America, thinning from above to low growing stock levels has even been found to result in a dramatic short-term increase in tree mortality (Powers et al. 2010).
Machine traffic causes damage to the roots of the remaining trees by exposing them to fungal attacks, resulting in financial losses due to discolouration and wood decay (Vasiliauskas 2001).On the other hand, the rutting of soil by harvesting machinery can positively affect CCF, as it exposes the soil, thereby promoting the emergence of new seedlings.According to Surakka and Sirén (2007) investigating the interlinkage between harvesting technique and silvicultural outcome in CCF is among the key topics requiring further research.Research into the seedling damage caused by the mechanised harvesting of boreal conifer forests is scarce and focuses on Norway spruce.Seedlings' height, their location relative to the remaining trees and strip roads, and harvesting intensity have been shown to significantly affect the probability of seedling damage in the selection cutting of Norway-spruce-dominated stands in Southeastern Norway and Eastern Finland (Granhus and Fjeld 2001;Surakka et al. 2011).
In the present study, we investigated the regeneration potential of Scots-pine-dominated mineral soil stands in CCF and RFM in Finnish Lapland, focusing on seedling damage.More specifically, the study aimed to quantify: (1) the total number of living seedlings that survived harvesting under different harvesting regimes (2) how different cutting methods and site preparation affect the proportion of damaged and undamaged seedlings The CCF cutting treatments included selection cutting and cap cutting; seed tree cutting and clearcutting were reference methods representing RFM.

Silvicultural treatments
The experiment was conducted in seven experimental forests (randomised blocks) located in the Scots-pine-dominated sub-xeric heath forests typical of southern and central Finnish Lapland (Figure 1).The experimental stands were harvested with the CTL method in October 2020.The treatments included clearcutting (CC), seed tree cutting (STC), gap cutting (GC), and selection cutting (SC).Four rectangular sample areas of 2 ha (100 m × 200 m) were established in each block and randomly subjected to the treatments above.These treatment areas were cut following the silvicultural guidelines of Metsähallitus, which is responsible for the management of state forests in Finland (Metsähallitus 2023).In the GC areas, four circular gaps with a diameter of 60, 40, or 20 m (2 gaps) were clearcut, and the intermediate areas were left untreated.On the SC plots, 50-60% of the basal area (BA) of Scots pine was removed.The target for the postharvest BA in the SC areas was 8 m 2 ha −1 , but the outcome varied due to the differences in stand structure.Most of the largest trees were removed in these areas, but some were left to secure adequate seed crops in the future.Dense groups of smaller trees were thinned.In the CC areas all merchantable trees were removed, while smaller trees, saplings, and seedlings were left.In the STC plots c. 50 seed trees per hectare were left.In addition, c. 10 retention trees per ha were left in groups in all treatments except for the GC plots.These retention tree groups covered c. 400 m 2 in each treatment area of 2 ha.The CC and STC areas were treated with disc trenching in the spring of 2021, except for two STC areas with a thin humus layer, following common practices in such stands.In addition, the undergrowth hampering soil preparation and forest regeneration was cleared before soil preparation in the CC and STC areas in one block, as recommended by Metsähallitus (2023).In the other six blocks, clearing was unnecessary.In the clearcut areas, seeding was integrated with disc trenching.
Five blocks were harvested with a Ponsse Scorpion King harvester and Komatsu 840 TX forwarder.In the other two blocks, a John Deere 1270E harvester and a Ponsse Buffalo forwarder were used.All the machines were equipped with tracks.The standard mass of the harvesters varies from 19.8 to 23.2, and that of forwarders from 13.8 to 19.8 t (Lectura 2023; Ponsse Plc 2023).Three harvester operators were involved in the cutting experiment, and each treatment area within a block was cut by the same operator.The John Deere harvester was used by one operator, and the Ponsse harvester by two operators.All had previous experience of the cutting methods included in the experiment.The harvesting conditions were normal for the season: the temperature ranged a few degrees Celsius both sides of zero, and in the latter half of October, the ground was covered by a 10-15 cm snow layer.In soil scarification a disc trencher attached to a forwarder was used.When pulling through soil preparation sites, the rotating discs trenched two parallel furrows approximately 180-200 cm apart.The discs formed 60-80 cm wide berms by casting organic matter, mineral soil, and slash to the outer edges of the furrows.A furrow network of 4000-5000 m ha −1 is considered to create favourable conditions for the emergence of a sufficient number of seedlings (Metsänhoidon suositukset 2023).

Field inventories
The stand volumes (stem volumes over bark) used in the modelling are based on postharvest relascope measurements and harvester data.Basal areas including trees with a minimum breast-height diameter (DBH) of 4.5 cm were measured by tree species from 10 relascope sample plots in each 200 m × 100 m treatment area.For each species, a representative tree was selected for height measurement with a clinometer, and the stand volumes (m 3 ha −1 ) were obtained, as described by Ärölä (2018).The initial stand volumes were calculated by summing the postharvest stand volumes above and the harvested volumes of merchantable timber (removals) based on the .hprfiles from the harvesters' production reporting systems (Arlinger et al. 2020).The mean stem volumes of removed trees for each treatment area were calculated by dividing the merchantable stem volumes by the numbers of removed trees obtained from the .hprfiles.
Damage to the seedlings was inventoried in the summer of 2021 from 25 circular sample plots of 50 m 2 (radius 3.99 m) from 3-5 lines in the 100 m × 200 m treatment areas.The inventory lines were placed at a distance of 20-45 m from each other perpendicular to the strip roads, and the sample plots were established on these lines with an interval of 20 or 24 m, as illustrated in Figure 2. The variable distances between lines and sample plots were used to avoid any systematic error harvester or forwarded operations within a treatment area could cause.In addition, before the cuttings eight permanent circular sample plots of 50 m 2 each were established per treatment area to monitor postharvest stand development (Natural Resources Institute Finland 2023).In these permanent sample plots, the mean number of seedlings before harvesting was on average (mean ± S.E.) 2818 ± 128 per hectare for CF, 2454 ± 113 for GC, 2957 ± 174 for SC and 2957 ± 215 for STC.
In the inventory, seedlings with a minimum height of 10 cm and breast-height diameter (DBH) less than 4.5 cm were included in the counts.Hence, both seedlings and saplings were recorded, but in the forthcoming text, both are referred to as seedlings.Only mechanical damage to the seedlings was considered.The status of the seedlings was categorised as follows: (1) undamaged; (2) slightly damaged; and (3) seriously damaged or dead.In addition, the coverage (% of ground cover) of the double humus layer reversed in soil scarification, slash (harvesting residue, i.e. tops and branches of harvested trees), exposed mineral soil, and strip roads were assessed in the circular sample plots.

Binomial model for the proportion of damaged seedlings
A logistic mixed-effects repeated model was applied to explain (and predict) the probability of damaged seedlings.The hierarchy in the model was a circular sample plot within a treatment area within a block (referred to by the indices i, j, and k).When the probability was applied to the seedlings on the 50 m 2 circular sample plot, it could be considered to be the proportion of damaged seedlings on a plot.Thus, the proportion was conceived of as the outcome of multiple binomial trials.In the binomial model, the number of damaged seedlings was treated as the number of successes, and the total number of seedlings minus the number of damaged seedlings as the number of failures.The data consisted of the counts of all the seedlings and damaged seedlings, and the multiple binomial trial was written as follows: ln number of successes number of failures = ln number of damaged seedlings total number of seedlings − number of damaged seedlings (1) The binomial model using random block effects was written as follows: where π is the probability of the event, i.e. the proportion of damaged seedlings.Binomial (n,p) denotes the binomial distribution with parameters, n describing binomial sample size, and p the proportion (probability) of the event's occurrence, i.e. seedling damage.
is a logit-link function, and f(.) describes the linear function with arguments X ijk (i.e.fixed predictors) and b (i.e.fixed parameters).The term m i describes the block variance, and m j describes the variance of the treatment areas.
Level k describes the circular sample plot level.A dispersion parameter in the logistic model was estimated in the model.The R package MASS and its function glmmPQL was used in the modelling (Venables and Ripley 2002).The (pseudo) R 2 values were computed using the R package MuMIn (Bartoń 2020).The R 2 values were computed for the marginal models (fixed predictors only used in computation).
Negative binomial models for the seedling counts Models 3-5 were computed to describe the total postharvest number of observed seedlings in the 50 m 2 circular plots and the postharvest number of undamaged seedlings.The latter represents undergrowth that could be utilised as complementary seedlings in forest regeneration.
A generalised linear mixed model with a negative binomial distribution assumption and log-link function was constructed to model the numbers of seedlings on the 50 m 2 circular sample plots (later referred to as seedling density models).The models computed using the R package glmmTMB (Brooks et al. 2017) consisted of three hierarchical levels: block, treatment area and circular sample plot, similar to the previously described binomial model.
In the R function, glmmTMB the variance function was defined in two ways: NB1 variance = m(1 + a) and NB2 Brooks et al. 2017).The NB1 parametrisation suggested a linear mean-variance relationship, and the NB2 parametrisation indicated a quadratic relationship.These two parametrisations were tried, and the better parametrisation (NB2) was selected based on Akaike's Information Criterion (AIC) and the models' simulated fits.The simulated distribution based on the model parameters theta (also referred to as size) and expected value (mu) were plotted against the observed count distribution to determine the overall model fit.
Generalised linear mixed models tend to underestimate the mean.A ratio estimator was therefore computed, and the estimated mean in the original scale was corrected with the estimated coefficient in the binomial and negative binomial models (Snowdon 1991).
The R package glmmTMB also allowed zero-inflation modelling.The possible existence of a zero-inflation problem was tested in the count models, but it was not an obvious problem in the data, based on the test with the R package DHRMa (Hartig 2021).The test with the testZeroInflation function indicated that there were no statistically significant differences between the observed and simulated zeroes (p ≥ 0.005).The (pseudo) R 2 values for the model were computed using the R package performance (Lüdecke et al. 2021).

Explanatory variables and model stability
The tested potential explanatory variables for the proportion and count models for seedling damage were the initial stand volume (stem volume, m 3 ha −1 ), removal (m 3 ha −1 ), the number of removed stems (pcs.ha −1 ), the average stem volume of removed trees (m 3 ), and the coverage (cover percentage, %) of slash, double humus, exposed mineral soil, and the strip road on the circular sample plots.Neither the postharvest stand characteristics nor the results from the seedling and soil coverage inventories were used in the same models, containing the cutting method as an explanatory variable due to the high influence of the cutting method on the postharvest stand characteristics.Instead, the initial stand characteristics could be tested in the same model with the cutting method.All the possible two-way interactions between the candidate variables were also tested and included in the models if they were statistically significant.
The differences in the mean values of the potential explanatory variables at block or circular sample plot level (Figures 3 and 4 in Data description) by the cutting method were tested using log-normal linear models with the random effects included in the models similarly to the proportion and count models.These simple models were computed to achieve unbiased standard errors and p-values for the differences between the cutting methods.The predictions of the models with their 95% confidence intervals were computed using the R package ggeffect (Lüdecke 2018) and the prediction graphics drawn using ggplot2 (Wickham 2016).
The risk of multicollinearity was considered in the model construction by checking the model stability (coefficients and p-values) after adding and removing the candidate variables for the models and checking the variance inflation factor (VIF) with the R package performance (Lüdecke et al. 2021).The VIF values of the main effects had to be less than five to be accepted in the same model as explanatory variables.

Data description
The potential explanatory data for the modelling were gathered at two hierarchical levels, as previously described.Preand postharvest stand volume, the number of removed stems and the average stem volume of removed trees were defined for the treatment area (2 ha) level, while the proportions of double humus, exposed mineral soil, slash, and the strip road were recorded for the 50 m 2 circular sample plots.
All the means of the variables illustrated in Figures 3 and 4 differed significantly by cutting treatment (p < 0.001).The preharvest stand volume was greater in the clearcutting and selection cutting areas than on the gap and seed tree cutting sites (Figure 3(a)).The original stock volume was therefore tested in the models as a covariate and included if it was significant at a 5% risk level.
The removal (m 3 ha −1 ) and number of removed trees per ha were greatest in the clearcutting stands, and the smallest in the gap cutting stands (Figure 3(b,c)).The proportion of the harvested area was greatest on the clearcutting sites and smallest on the gap cutting sites, where the gaps constituted only c. 40% of the total treatment area.This should be kept in mind when interpreting the modelling results.The largest trees were harvested from the selection cutting areas, and the smallest from the gap cutting areas (Figure 3(d)).Logging residues (slash, Figure 4(a)) and strip roads (Figure 4(d)) covered more than 20% of the area on many circular sample plots.The mean and median values of the summed soil disturbance factors (proportions of exposed mineral soil, double humus layer, slash and strip roads) on the circular sample plots were 15% and 16% respectively, ranging mostly from 0% to 40%.The total soil disturbance differed from 0%, especially in the clearcutting and seed tree stands.The raw data used in these computations suggested that most of the circular plot area remained undisturbed.

The proportion of damaged seedlings
The proportion of damaged seedlings varied markedly between cutting treatments (Figure 5, Model 1 in Table 1).The greatest proportions of damaged seedlings were encountered in clearcutting and seed tree cutting, where 25-34% of the seedlings were damaged (Figure 5(a)).The predicted proportion of damaged seedlings in gap cuttings was c. 7% and 12% in selection cutting.These estimates were computed using the average preharvest stand volume (112 m 3 ha −1 , Figure 5(b)).
The proportion of damaged seedlings increased with increasing preharvest stand volume (Model 1 in Table 1, Figure 5(b)).This proportion increased from 23% to 66% when the preharvest stand volume increased from 77 to 191 m 3 ha −1 (Figure 5(b)).These predictions were computed by fixing the other predictors to the level of "clearcutting".
The cutting method was excluded from the other binomial model (Model 2 in Table 1, Figure 6) for the proportion of damaged seedlings.Instead, the continuous covariates at forest stand or sample plot level were included.The model shows that the coverage of strip roads had the strongest effect on the damage proportion of seedlings (Table 1, Figure 6(a)).If half the sample plot area was covered by strip roads, about 90% of the seedlings were damaged, and all the seedlings were damaged when the coverage of the strip road exceeded 75% (Figure 6(a)).When the coverage of exposed mineral soil, double humus or slash (logging residues) increased from 0% to about 25%, the proportion of damaged seedlings increased from 10% to 30-60% ).In addition, an increase in removal at the treatment area level slightly increased the occurrence of seedling damage.When the removal increased from 25 to 180 m 3 ha −1 , the proportion of damaged seedlings increased from c. 10% to 25%.

The number of undamaged seedlings
There were only minor differences between the cutting methods in the total postharvest numbers of seedlings (Model 3 in Table 2, Figure 7(a)) computed using the mean value of the preharvest stand volume (128 m 3 ha −1 ).However, significant differences in the numbers of undamaged seedlings were found (Model 4 in Table 2, Figure 7 (b)).There were about 1000 undamaged seedlings ha −1 in the gap cutting and seed tree areas after harvesting, and about 750 undamaged seedlings ha −1 in the clearcutting areas.The number of undamaged seedlings in the selection cutting areas was twice that of the clearcutting areas (Figure 7(b)).Higher preharvest stand volume significantly reduced the number of undamaged seedlings after cutting (Model 4 in Table 2, Figure 7(d)).
The cutting method was replaced by the statistically significant continuous variables (covariates) in the second version of the negative binomial count model for the undamaged seedlings after cutting (Model 5 in Table 2).The results show that none of the measured forest characteristics (removal, the average stem volume of removal, the preharvest number of stems), or the coverage of exposed mineral soil significantly affected the postharvest number of undamaged seedlings (p > 0.05).On the other hand, increasing the coverage of the strip roads, double humus, and slash considerably decreased the number of undamaged seedlings on the circular sample plots (Figure 8(a-c)).In particular, strip roads had a detrimental effect on the number of undamaged seedlings.As the coverage of the strip roads approached Table 1.The generalised linear mixed effects models for the proportion of damaged seedlings with a binomial distribution assumption.Std.error denotes the standard error of the estimate, df denotes degrees of freedom, t-value denotes the t-value for the parameter estimate, chi-value denotes the chi-squared test value for the Type III Anova test for the categorical variable, and CI denotes confidence interval.100%, the number of undamaged seedlings approached zero (Figure 8(a)).

The performance and fit of the models
In general, the (pseudo) R 2 values of the binomial and negative binomial models remained rather low (from about 16% to 20%), the highest R 2 value reaching 35% with the five soil disturbance factors.
The binomial models for the proportions of damaged seedlings fitted moderately to the data when the distributions of the observed and predicted values were compared.Again, the model with five continuous predictors fitted considerably better than the model including only the cutting method and the postharvest stand volume (Table 3).
The negative binomial models underestimated the mean value of the postharvest number of seedlings, similarly to the binomial models.However, the simulated distributions of the negative binomial models compared to the observed distributions (Appendix) suggested a fair fit of these models to the data.The negative binomial count models for the number of undamaged seedlings estimated the count distribution from 1 to about 20 seedlings on a sample plot (50 m 2 ) quite well but underestimated the number of zeroes and the highest counts.The model for the total postharvest number of seedlings was clearly best fitted to the data.

Discussion
In the present study, the greatest proportions of damaged seedlings (25-34%) were found in the areas harvested with the rotation clearcutting and seed tree methods.In the CCF methods, the damage proportions were much lower: 7% in gap cutting and 12% in selection cutting.However, the cutting treatment as such was not a very strong predictor of the seedling damage.Harvesting machinery traffic and soil scarification in the case of clearcutting and seed tree cutting were the most important causal agents of the seedling damage.Accumulated slash on the seedlings also damaged them.The results are in accordance with the study of Hyppönen (2000), who found that an increase in removal and strip road coverage increased the proportion of damaged Scots pine seedlings during the cutting of hold-over trees in Lapland, although the preharvest stand volume was more important contributor to seedling damage than removal in the present study.
The initial number of seedlings in the sample plots established for assessing seeding damage are not known.However, based on the preharvest seedling counts (2400-2957 pcs ha −1 ) in the permanent sample plots for monitoring stand development (Natural Resources Institute Finland 2023), it can be concluded that significant amounts of seedlings were lost due to soil disturbance and harvesting operations.In this study, binomial modelling was applied, i.e. the seedlings were considered to be either undamaged or   damaged (0/1), although the severity of the damage was assessed in the field inventory.In practice, some of the damaged seedlings still had development potential.The type of seedling damage or the size of seedlings was not recorded.In the study of Surakka et al. (2011), stem breakage was the most common type of injury for 0.5-2.5 m high spruce seedlings on sites harvested with selection cutting.
In RFM, the seedling damage caused by harvesting and soil scarification is not detrimental, as clearcut areas are either sown or planted, and in seed tree cutting, natural regeneration is expected after cutting and site preparation.Moreover, soil exposure by harvesting machinery usually promotes the emergence of natural seedlings.Due to its assumed low technical quality and growth potential, undergrowth is often cleared after harvesting to avoid an uneven stand structure.Hence, in periodic forestry based on clearcutting or seed tree cutting with subsequent site preparation, there has been no reason to save natural seedlings and other undergrowth from the forest regeneration perspective.On the contrary to RFM, the number of seedlings is an important factor to consider in stand management using CCF methods, particularly in northern areas where a cool climate restricts natural seed production and, consequently, stand regeneration (Juntunen and Neuvonen 2006;Hilli et al. 2008;Pukkala et al. 2010).
The estimates for preharvest stand volumes, obtained by summing the postharvest stand volumes and the removals based on the harvester data, are underestimations.This is because the data on harvested production did not include unmerchantable stem volumes, which consist of undersized tops, reject sections, and cleared undersized stems.The proportion of unmerchantable top sections depends on tree size: the smaller the tree is, the greater the proportion of unmerchantable stem volume (Hakkila et al. 1995).However, the bias described above probably does not affect the conclusions drawn in the present study, as the error is systematic, i.e. similar in all treatments, and the experimental stands were relatively homogenous.
CCF is based on the presence, survival, and growth of understorey seedlings and saplings, and trees of lower canopy layers (Surakka et al. 2011).Size and competition status are significant factors affecting the future development of undergrowth (Eerikäinen et al. 2007).Especially in selection cuttings, the competition with the remaining upper storey trees can inhibit seedling establishment and survival.Avoiding damaging to the seedlings with forest machinery is therefore very important.In this study, the average number of undamaged seedlings in the selection cutting areas (ca.1500 pcs ha −1 ) was satisfactory, but their spatial distribution was not recorded.The number of viable seedlings may therefore be overestimated, as some of the seedlings will be suppressed in the future (Kyrö et al. 2022;Rautio et al. 2023).In uneven-aged stands, seedlings often grow in groups, and there are open areas lacking seedlings.Hence, from the spatial distribution perspective, damaging a single seedling can be more harmful than losing some of the clustered seedlings (Surakka and Sirén 2007).On the other hand, complete seedling coverage is not required in CCF as in the case of clearcuts in RFM.However, it has been shown that in stands with a high basal area, the average seedling growth is extremely low, indicating that older seedlings are dying, and new seedlings are emerging all the time (Niemistö et al. 1993;Rautio et al. 2023).Consequently, ample seedling material recorded in a given inventory may merely change constantly instead of developing into more mature trees.
Avoiding damage to the remaining trees is also important for reducing the risk of root rot (Heterobasidion parviporum or H. annosum) infections in CCF (Hantula et al. 2022).In the study of Eliasson et al. (2003), the probability of sapling breakage in harvesting increased with sapling height and decreasing temperature.They recommended that harvesting during a severe frost (< − 15 °C) should be avoided as the risk of damage to the remaining stand is high (Eliasson et al. 2003).On the other hand, frost prevents soil rutting and thereby reduces damage to the roots of the remaining stand, controlling the risk of fungal attack (Shoop 1995;Saarilahti 2002).Furthermore, snow protects small seedlings.Targeting the operations at winter months with a mild temperature could therefore be a good strategy for avoiding seedling damage.Additionally, Heterobasidium sp.does not produce and disperse spores at low temperatures, which makes the cold season a safe harvesting time.However, the timing of harvesting may become operationally challenging if stand management using CCF methods increases markedly.
On clearcut and seed tree sites, the purpose of soil preparation is to promote stand establishment.The aim is therefore to treat the area systematically.This increases the risk of seedling damage if the information about the seedling clusters to be saved does not reach the person operating the soil preparation machine.Modern harvesters collect a lot of georeferenced data.To protect the seedlings or saplings to be saved, their location could be recorded while cutting.In RFM, the use of harvester data has been suggested in precise tree species selection within stands infected by root rot (Saksa et al. 2021).The cutting of trees did not seem to cause much damage to seedlings in the CCF stands, indicating that cutting techniques can be adjusted to CCF with reasonable silvicultural outcome.Besides the productivity of harvesting, the operator has a crucial impact on the silvicultural outcome of harvesting.Consequently, the differences in operators' skills are emphasised in the harvesting of uneven-aged stands (Surakka and Sirén 2007).The coverage of strip roads and the risk of seedling damage can be minimised by careful planning of strip roads.A stable base machine with high ground clearance and a long-reach crane enables wide strip road spacing in cutting and forwarding.Sufficient capacity to handle large trees above seedlings is also required of the machines used in selection cutting.Decreasing the width of tyres or tracks of the machines to reduce seedling damage is not a realistic option for reducing seedling damage if the forest machines remain as large as they are now.

Figure 1 .
Figure1.Location of the experimental blocks (n = 7) and an example of the experimental design with four 2 ha treatment areas randomly subjected to the cutting treatments.

Figure 2 .
Figure 2. The location of circular sample plots established for the inventory of seedling damage and soil disturbance within the treatment areas.

Figure 3 .
Figure 3. Boxplots describing the distributions of the stand characteristics by cutting method in the 2 ha treatment areas: the initial stand volume (m 3 ha −1 ; a); removal (m 3 ha −1 ; b); number of removed trees (pcs.ha −1 ; c); and the average stem volume of a removed tree (m 3 ; d).

Figure 4 .
Figure 4. Boxplots describing the coverage (% of ground cover) of slash (a), the double humus layer reversed in the soil scarification (b), exposed mineral soil (c), and the strip road (d) on the 50 m 2 circular sample plots.

Figure 5 .
Figure 5.The predictions (with 95% confidence intervals) of the binomial model for the damaged proportion of seedlings by cutting method (a) and the preharvest stand volume (b).The predictions for the cutting methods were computed for the mean preharvest stand volume (112 m 3 ha −1 ), and the predictions for the preharvest stand volume at the level of "clearcutting".

Figure 6 .
Figure 6.The predictions (with 95% confidence intervals) of the binomial model for the damaged proportion of seedlings by the coverage of the strip road (a), exposed mineral soil (b), double humus (c), slash (d), and removal (e).The values of the other variables were set to their means.

Figure 7 .
Figure 7.The predictions with their 95% confidence intervals of the negative binomial model for the total postharvest number of seedlings (a) and the postharvest number of undamaged seedlings (b) by cutting methods and the preharvest stand volume (c-d).The predictions for the cutting methods were computed at the mean value of the preharvest stand volume (128 m 3 ha −1 ), and the predictions for the preharvest stand volume at the level of "clearcutting" in the variable cutting method.

Figure 8 .
Figure 8.The predictions (with 95% confidence intervals) of the negative binomial model for the postharvest number of undamaged seedlings by the coverage (% of ground cover) of the strip road (a), double humus (b), and slash (c).The values of the other variables (than those predicted) were set to the mean.

Table 2 .
The generalised linear mixed effects models for the number of undamaged seedlings with a negative binomial distribution assumption (quadratic variance function parametrisation).Std.error denotes the standard error of the estimate, df denotes degrees of freedom, t-value denotes the t-value for the parameter estimate, chi-value the chi-squared test value for the Type III Anova test for the categorical variable, (the degrees of freedom in parenthesis), and CI denotes confidence interval.

Table 3 .
Observed and predicted proportion distributions of the binomial models for the damaged proportion of seedlings.