Forests in the study were managed by Military Forests and Estates of the Slovak Republic, Forest district Kežmarok (Figure 1; GPS coordinates: 49°15′12.5″N 20°31′47.98″E). The Forest district was located in the northern part of Slovakia, in the steep, rugged slopes of Levoča mountain range (Western Carpathians).

The Forest district managed 17,833 ha of forests with elevation between 600 and 1,290 m asl. The most abundant species was Picea abies L. (63.2%), followed by Fagus sylvatica L. (14.4%), Pinus sylvestris L. (9.1%), Abies alba Mill. (5.8%), Larix decidua Mill. (3.5%), with annual harvests reaching ~45,000 m³ of timber and a forest road network density of 8.21 m/ha.

Due to the size and diversity of the area (from the perspective of the natural conditions), we only report on soils in stands with the measurement plots (Table 1). Umbric Podsols Schad (2016) created from lighter solid products of weathering of acidic bedrock were the soil type in observed forest stands. According to soil particle size and texture classes, the soils were sandy loams, with a moderate soil water retention capacity (95–130 l/m²) and a moderate soil infiltration capacity (0.081–0.8 cm/min). The bedrock was created mainly from flyshoid sandstones. Mean annual rainfall in the observed mountain region ranges between 800 and 900 mm, and the mean annual temperature varies between 4 and 6°C. The mean potential evapotranspiration ranges between 400 and 450 mm. Detailed climatic data from the year 2022 were recorded by a GARNI 1025 Arcus (GARNI technology, Ostrava, Czech Republic) weather station located in Lubica village (GPS: 49°07′08.9″N and 20°29′42.9″E), the closest (16 km) weather station to the study area. According to the weather station, 595 mm of rainfall (monthly: Ø 49.55; min. 11.40; max. 107.72) occurred in 2022, which is slightly above normal (12 mm) compared to long-term rainfall (583 mm) at the given weather station. However, during August and September, when the harvesting and field work were conducted, the area experienced a rainfall deficit (August: −25.25%; September: +5.02%).

The machinery was used as follows: a ROTTNE H11^C harvester and a ROTTNE F15^D forwarder were used to fell, limb, buck, and extract timber from the forest stand to the roadside while secured by two Ecoforst TWINCH 10.1 mobile traction winches (i.e., a single winch per machine), located on the hilltops and radio-controlled by the harvester and forwarder operators (Table 2). The harvester carried out thinning operations and passed over 4 m wide trails, which the operator created (Supplementary Figure 1). Winch-assisted harvesting took place in stands with steep slopes (46.4–51.8%) no. 1, 2, and 3. The CTL machinery worked without a winch-assist in the Control stand because the slope was not steep (26.7%) and enabled safe traffic (Table 1).

We conducted the measurements between October 10 and 15, 2022, immediately after the harvesting machinery left the forest stands, and it was safe to begin measuring. During fieldwork, no precipitation occurred, and temperature ranged between 7 and 15°C. Forest stands were selected based on their similarity in soils, slopes, tree species mix, tree spatial distribution, other natural conditions, and proximity (Table 1; Figure 2). Sample plots were established on particular forwarding trails. The sample plot size was 10% of the harvested area in stands up to 50,000 m² or 5% in stands where the harvested area was larger than 50,000 m² (Allman et al., 2015). Plots were squares, with a side of 20 × 20 m. The positioning of the sample plots was calculated according to Scheer and Sedmák (2010) (Equation 1):

s, sample plot distance (m)

P, harvested area (ha)

n, number of sample plots

Measurement locations were established on the opposite sides of each plot (Figure 3). In each measurement location, we (i) determined the longitudinal slope with a True Pulse 360B (Laser Tech Inc., Colorado, USA) laser range finder (in °); (ii) evaluated the presence and depth of slash on the trail (cm), (iii) measured the forwarding trail profile (cm) and maximal rut depth, (iv) collected soil samples from the localities stand, rut and in the middle of the trail to determine the bulk density (BD) of soil (g cm⁻³) and soil moisture content (%), (v) measured soil penetration resistance (PR) in MPa and penetration depth (cm).

To measure the average forwarding trail profile depth (cm), we placed a 4 m ranging rod perpendicularly on the trail. We measured the trail profile with a tape measure from the soil surface to the lower side of the rod, each 20 cm of the length of the ranging rod (Jankovský et al., 2019). Maximal rut depth (cm) was expressed as the largest value measured in the ruts of the forwarding trail. The intensity of disturbance via rut creation was classified according to the rut disturbance classification scale, shown in Table 3. Soil samples were collected into 100 cm³ Eijkelkamp sampling cylinders (Giesbeek, Netherlands). Immediately after sampling, the cylinders were sealed in a plastic bag for subsequent BD and moisture content laboratory analyses. BD of soil and moisture content (Equation 2) was determined in a dried state from samples dried for 24 h at 105°C (Allman et al., 2015).

Penetration resistance and penetration depth of soil were measured in the exact locations of stand, rut and middle by an Eijkelkamp Penetrologger (Giesbeek, Netherlands), equipped with an 80 cm extension rod and a 30° cone screwed onto the bottom end of the rod. The cone had a nominal diameter of 11.28 mm and a base area of 1 cm². Thirty-nine sample plots containing 78 measurement locations were established in total. From the plots, 234 BD and moisture content samples were collected, and 234 PR measurements were carried out. Penetration resistance and penetration depth were measured by an Eijkelkamp Penetrologger, which recorded values per 1 cm of penetration in the range between 0 and 80 cm, based on the depth to which we were able to push the penetrometer which resulted in a higher number of values for these variables.

Data were processed in Tibco Statistica 14.0.1 and IBM SPSS statistical software. Shapiro-Wilk and Levene's tests were used to test the normality of data distribution and the homogeneity of variance for BD, moisture content, PR, and penetration depth.

One-way analysis of variance (ANOVA) was used to compare the differences between the forest stands (Duncan's test) and measurement locations (Tukey's test) for BD, soil moisture, PR and penetration depth, as well as to compare the differences between the mean and maximal rut depths formed in the observed forest stands. Linear mixed effects models (LMMs) were used to identify the effects of a particular forest stand (random effect) and measurement location (fixed effect) on the variability of BD, soil moisture, PR, and penetration depth, as well as to identify the effects of the presence of traction winches, slash in the forwarding trails, slope, and soil moisture content as fixed effects on BD, PR, and maximal rut depth in forwarding trails as dependent variables. Furthermore, we used LMMs to observe the relationships between PR and BD and the thickness of the slash cover and PR and BD on data from the rut location. We hypothesize that CTL machines cause significantly higher H₁ soil compaction in the ruts, expressed as BD; H₂ penetration resistances (PR) in the ruts of the forwarding trails; H₃ rut depths in forwarding trails (cm), compared to CTL machines deployed in standard conditions without winch assistance. The hypotheses H₁, H₂, and H₃ were tested on a 5% α level of significance, i.e., they were rejected, and the alternative hypothesis H_A “no significant differences between the variables were confirmed” was accepted.

From the mean dried BD values, it can be seen that machine traffic increased BD in the rut by 13.42 (no. 2)−22.35% (Control), representing the most substantial BD increase from all observed forest stands. In contrast, BD was lower in the middle (0.77–0.86 g.cm⁻³) than in the stand (0.82–0.90 g.cm⁻³), probably due to the soil displacement into the trail's center caused by machine traffic. The lowest BD values in the ruts were recorded in stands no. 2 (0.93 g.cm⁻³) and Control (1.04 g.cm⁻³), where the lowest amount of timber was harvested (221.79 and 109.34 m³, respectively), resulting in smaller passages of forest operations over the trails. The highest values of moisture content (%) were recorded at the rut location in all forest stands, except for stand no. 3, and fluctuated in the range of 44.24 (no. 1)−54.29 (no. 2) % (Table 4). The machine passage and the resulting change in the microrelief in the trails caused higher moisture content (%) in the ruts. One-way ANOVA (Duncan's test) did not confirm significant differences in BD between the forest stands (no. 1, 2, 3, Control) at the exact measurement locations, i.e., stand (F = 0.638; p = 0.593), rut (F = 2.373; p = 0.077), and middle (F = 1.553; p = 0.208). The results showed that the conditions and the intensity of soil disturbance in the ruts and middle locations were similar between the forest stands. The differences between locations (Tukey's test) were not significant only in the case of the Control stand (F = 2.587; p = 0.09) (Table 4). Regarding moisture content, the differences between the four observed forest stands were insignificant only at the rut location (F = 2.503; p = 0.07). Moreover, in the case of stand no. 3, the differences between the stand, rut, and middle measurement locations were insignificant (F = 14.816; p = 0.00) (Table 4).

Linear mixed model of the effects of measurement locations (stand, rut, middle) as fixed effects, and the forest stand no. (random effect) on BD showed that 18.9% (Rm2 = 0.189) of the BD variability was affected by measurement location (marginal effect) and 22.4% (R²c = 0.224) by fixed and random effects (conditional effect). Only about 3.5% of the variability could be attributed to random effects. Measurement location as the fixed effects significantly affected BD (F = 28.297; p < 0.05). In the case of moisture content as a dependent variable, the R²m was 0.067, and R²c was 0.202, and the random effects (stand no.) explained more variability of the moisture content (13.5%) than fixed effects (6.7%) and had (fixed effect) a significant relationship with moisture content (F = 9.749; p < 0.05).

We then observed the relationship between the BD in ruts and the usage of traction winches, slope category, and moisture content category as fixed effects and the stand no. as random effects via an LMM (Table 5). The analysis showed that 73.2% (Rm2 = 0.732) of rut BD variability was explained by fixed effects and 75.7% (R²c = 0.757) by the combination of fixed and random effects (Table 5). The random effect (stand no.) contributed only 2.5% to explaining rut BD variability, although only moisture content significantly affected the dependent variable (F = 49.164; p < 0.05).

Penetration resistance increased in the rut location, from 4.88 (no. 1) to 26.8% (Control). Table 6 shows that the increase in penetration resistances between the rut and stand locations was largest at the Control stand. Compared to the stand location (2.50–3.10 MPa), the PR in the middle of the trails (2.57–3.06 MPa) varied in some cases, it increased by 8 (Control) to 17 % (no. 2), while in other cases, it decreased by 1.29 (no. 3)−10.45 % (no. 1). There was a decrease of penetration depth (cm) 1.96 (no. 3)−10.95% (no. 2) in ruts (18.98–19.77 cm) compared to the stand location (19.36–22.20 cm) because the cone's ability to penetrate the soil was decreased as a result of machinery traffic. Interestingly, in stand no. 1, the penetration depth (cm) increased by 33.48% (4.8 cm) in the ruts compared to the stand location (Table 6). At the middle location, penetration depth (cm) increased by 0.26 (no. 3)−17.83% (no. 1), compared to the stand location, except for stand no. 2, where we observed a decrease in penetration depth by 7.79%. One-way ANOVA (Duncan's test) confirmed that the differences between the forest stands (no. 1, 2, 3, Control) were significant on the stand (F = 19.29; p = 0.00), rut (F = 7.759; p = 0.00), and middle (F= 54.367; p = 0.00) measurement locations. The differences between the stand, rut, and middle measurement locations (Tukey's test) were also significant in all forest stands (p < 0.05). In the case of penetration depth (cm), Duncan's test confirmed that for the same measurement locations, the differences between forest stands (no. 1, 2, 3, Control) were significant (p < 0.05), except for the rut location (F = 1.368; p = 0.25). Tukey's test showed that the differences between the stand, rut, and middle measurement locations were significant in all stands except for stand no. 3 (F = 0.189; p = 0.828) (Table 6).

Linear mixed modeling between PR as the dependent variable, measurement location as fixed effects, and stand no. as random effects showed that only 0.7% (Rm2 = 0.007) of PR variability were affected by measurement location, and further 4.3% were affected by fixed and random effects (R²c = 0.05) and the relationship was insignificant (F = 0.446; p = 0.64). We achieved similar results when penetration depth (cm) was the dependent variable (R²m = 0.005, R²c = 0.012). The LMM between the rut PR, stand no. (random effect), and the usage of traction winch (F = 2.222, p = 0.164), presence of slash in the ruts (F = 3.134, p = 0.081), slope category (F = 0.913, p = 0.464), and moisture content category (F = 1.704, p = 0.160) (fixed effects) showed no significant relationship (Table 5). Another LMM did not confirm a significant relationship between PR and BD (F = 0.901; p = 0.346; R²m = 0.012).

Slash (treetops and limbs) were placed on the trails by the operator during the harvesting process to minimize the effects of the machine on the soil surface. From measured data, it can be seen that the share of measurements without slash placed on the trails was relatively high, ranging from 67 to 95% (Table 7). Consequently, LMM in the ruts of the forwarding trails did not confirm a significant relationship (F = 0.264; p = 0.769) between the slash thickness category (cm) and BD (R²m = 0.007; R²c = 0.094). In the case of PR and slash thickness category (cm), the relationship was significant (F = 3.819; p = 0.027), though the explanatory power was weak R²m (0.090) and R²c (0.094). Insignificant or weak relationships between slash cover and BD or PR were caused by the small amount of slash cover over the ruts of the forwarding trails, which ranged between 4.17 and 29.16%. An inadequate layer of slash on the trail or its absence was mainly caused by the harvesting operations being carried out in young, 30- to 40-year-old stands. The thin branches and the low harvesting concentrations did not permit covering the trails with sufficient slash.

Mean forwarding trail profile depth in particular forest stands varied between 4.63 and 7.28 cm, whereas maximum rut depth ranged between 10.86 and 17.25 cm (Table 8). The deepest average trail profile and maximum rut depth were recorded in the Control stand.

The most frequent intensity class of forwarding trail disturbance was a moderately strong disturbance, with a maximum rut depth between 8 and 15.9 cm (Table 8). This level of disturbance represented 77.27 (stand no. 3) to 80% (stand no. 2) of the trails in the stands. In the Control stand, severe disturbance was the most frequent, representing 50 % of disturbances, with a maximum rut depth between 16 and 25.9 cm. One-way ANOVA (Duncan's test) confirmed significant differences between stands. In the case of average rut depth (cm) (F = 7.21; p = 0.00) and maximum rut depth (F= 7.48; p = 0.00), the differences were significant in the Control stand, which differed significantly from all other forest stands (nos. 1, 2, 3).

The LMM for maximum rut depth (dependent variable), traction winch usage, slash cover, slope category, and moisture content as fixed effects, and stand no. as a random effect showed a significant relationship between the traction winch usage (F = 9.091; p = 0.019) (Table 5). In the case of other variables, i.e., slash cover (F = 1.231; p = 0.271), slope category (F = 1.724; p = 0.156), and moisture content category (F = 0.379; p = 0.823), the relationships were insignificant.

Garren et al. (2019) studied the changes in soil parameters when mobile traction winches are used in steep slopes (27–38°, i.e., 51–78%), dividing the slopes into three categories (first: 27–31°, i.e., 51–60%; second: 31–35°, i.e., 60–70%; third: 35–38°, i.e., 70–78%). To compare, we chose the first slope category from the study by Garren et al. (2019). In that category, the bulk density in rut ranged from 1.4 to 1.58 g cm⁻³ and increased by 16% on average compared to stand locations. In our case, at 22–29.5° (40.4–56.6%) slopes, the bulk density in rut ranged between 0.93 and 1.05 g cm⁻³ and increased by 13.42 (stand no. 2)−22.09% (stand no. 3). At lower slopes (14.9° i.e., 26.6% stand Control), higher bulk density was not observed in the rut (1.04 g cm⁻³) compared to other stands, though the increase was more substantial (22.35%) than in stands with higher slopes. The LMM did not confirm a statistically significant effect of slope on the BD (p > 0.05). To our knowledge, very few studies focus on the problem of using mobile traction winches combined with CTL machinery in thinnings and the soil disturbance caused by this harvesting system. However, we can use studies of non-assisted CTL machinery for comparisons. Malvar et al. (2017) indicate that mean BD at 0–5 cm depth ranges from 0.77 kg cm⁻³ at the stand area to 1.07 kg cm⁻³ at compacted rut area, and BD is between 14 and 28% higher in rut. This research was conducted on (47–53%, i.e., approximately 25.2–27.9°) slope. Allman et al. (2015) report an increase in soil bulk density from (0.82 to 1.12 g cm⁻³) in stand to (1.02 to 1.36 g cm⁻³) in rut, i.e., between 21 and 2% after CTL machine traffic. Labelle et al. (2019) indicate a pre-operation soil BD from 0.4 to 1.0 g·cm⁻³ and from 0.7 to 1.2 g cm⁻³ after machine traffic, representing an increase between 20 and 75%. Williamson and Neilsen (2000) report a BD increase of 0.17 g.cm⁻³ after machine traffic.

Penetration resistance increased by 4.8–26.8% in the rut, while the penetration depth decreased by 2–10%. The most substantial increase in PR (26.8%) was observed in stand Control, where the winch did not operate. However, according to Tukey's test, the differences between the stand and rut locations were significant in all forest stands. Labelle et al. (2019) state that on forest or agricultural land, high soil PR may be problematic for root growth and values beyond 3.0 MPa are regarded as a threshold for interrupted root growth, thus marking severely disturbed soils. In our case, PR over 3.0 MPa was recorded even in undisturbed locations, although only in the case of stand no. 3. To compare, PR values in the rut ranged between 3.01 and 3.17 MPa. Labelle et al. (2019) concluded that while none of the pre-forwarding measurement points exceeded the 3.0 MPa threshold, 10.7% of post-forwarding measurement points showed PR above 3.0 MPa, with about 80% of such values occurring in the rut. We presume that the higher PR values in the stand locations in our case were due to insufficient rainfall during August and September, as well as moisture deficit in lower soil layers. On the other hand, the LMM did not confirm a statistically significant relationship between the soil moisture content on the PR (p > 0.05), likely due to our use of determining the moisture content from samples taken from upper soil layers via a gravimetric method. In contrast, the PR values were recorded at depths up to 80 cm, where different moisture conditions can occur, which might have caused an increased PR when trying to penetrate the cone penetrometer into deeper soil layers. Cudzik et al. (2017) report even higher PR values. Their study showed that PR after CTL machine traffic increased by 220–308% compared to control measurements at the stand location. The penetration resistances measured ranged from 2.8 to even 4 MPa.

In their study, Garren et al. (2019) indicate that the maximum rut depth (cm) caused by the winch technology was 9.5 cm and in 39% of observations, the rut depth was 0 cm or no visible ruts and overall rutting was negligible. In our case, the maximum rut depth ranged between 10.86 and 17.25 cm. In more than one-quarter to one-fifth of observations where winch assistance was used, the soil disturbance reached a “moderately strong” category, while in the stand without a winch, the rutting reached a “severe disturbance” level most frequently, though the differences in rut depth (cm) were primarily insignificant. Visser and Spinelli (2023) state that on steeper slopes (30–40%, i.e., 16.7–21.8°), low soil disturbance and rut deeper than 15 cm covers only 2% of the observed surface after being trafficked by a two-axle skidder with winch assistance. In our case, ruts deeper than 15 cm were measured frequently, between 9 (stand no. 3) and 50% (Control) of measurement locations. Although traction-assist technology significantly reduces wheel slippage, no mitigating effect of traction-assist technology on soil disturbance was observed. Cudzik et al. (2017) report a rut depth of up to 5 cm on 81% of trails in their study of CTL machines working in slopes between 8 and 12° (14–21.3%). Han et al. (2006) state that ruts are deeper for moist (medium and high moisture contents) than for drier soils. The high moisture content (%) contributes to the development of ruts that are significantly (13.8–35.3 cm) deeper than in the case of medium (8.2–12.0 cm) or low moisture contents (4.5–8.0 cm). In our case, the deepest ruts were recorded in the forest stand, where the CTL machines operated unassisted, as confirmed by the fact that the LMM showed a statistically significant effect of a traction winch being used (p < 0.05) on the maximum rut depth in the forwarding trails. Duncan's test confirmed significant differences between the mean and maximum rut depth between the Control stand and stands no. 1, 2, and 3.

Forest harvesting in mountainous regions of the Carpathian mountain range happens in various technological conditions connected to machinery and their use, soil types, terrain configuration, and tree species being harvested. Therefore, large-scale studies are needed to cover the array of machines and their use in highly variable natural conditions. Future research should expand the use cases of traction winches supporting ground-based wheeled machinery and focus on repeated observations with control plots that reflect the variability of natural conditions. Different approaches to observing the soil disturbance complementary to the methods used in this study can also be used in the future.

Admittedly, our study has a limited scale. We could only perform control measurements in one comparable forest stand, while more would be preferable. Our inability to do so was caused by the fact that no suitable forest stands were being harvested by the land owner where the slope would be sufficient on the one hand, and occupational safety and health could be secured on the other since the operation of unassisted machinery in slopes over 40% slope was forbidden by the internal regulations of the land owner.

Nevertheless, even small-scale studies contribute to filling in the considerable knowledge gap on the environmental effects of this technological variant and providing practitioners with adequate information when considering what technology and machinery to deploy in their forests.