Introduction

Protected areas such as national parks, nature reserves, and other wilderness areas play an important role in conserving species that are unlikely to survive in intensively managed areas and ecosystems (Branquart et al. 2008). In current times, protected areas also play a significant role in the mitigation and adaptation of biological entities to climate change and in providing several ecosystem services. Well-managed protected areas reduce rates of habitat loss, can provide carbon sequestration, and reduce deforestation. There is also strong evidence that protected areas maintain species population levels, including threatened species, better than other management approaches (Lubchenco et al. 2003; Postel and Thompson 2005; Scharlemann et al. 2010; Soares-Filho et al. 2010; Watson et al. 2014). Unfortunately, so far, this evidence does not extend to soils and soil communities, justifying the need for more targeted approaches (Zeiss et al. 2022).

Despite the importance of soil biodiversity for ecosystem functions, defined as physicochemical and biological processes that take place in a given ecosystem, and that maintain terrestrial life, and services, that are ecosystem functions that directly benefit human life (Kremen 2005), nature conservation policies are still mainly concerned with aboveground species (Guerra et al. 2021a, b). Soil biodiversity has been largely ignored, resulting in a lack of information about the conservation status of most soil organisms and about the effects that nature conservation policies have on soil systems (Ciobanu et al. 2019; Guerra et al. 2022; Zeiss et al. 2022). Currently, one of the main strategies for nature conservation is the establishment and maintenance of protected areas; and a significant amount of time, effort, and resources have been invested in their prioritization and management in the last decades (Elbakidze et al. 2013; Gaston et al. 2008). In Portugal, the Institute for Nature Conservation and Forests (ICNF) is the mandated Portuguese governmental body for nature and forest policies, including the management of Protected Areas and State-managed national, municipal, and communal forests of mainland Portugal (ICNF 2017). Forest areas are particularly important for conservation because of their use as a refuge and habitat for biodiversity and as ecosystem service suppliers (Cavard et al. 2011; Demissie et al. 2019; Leberger et al. 2020; Raina et al. 2011). They play an important part in the mitigation and adaptation to climate change as trees have important roles in carbon sequestration and storage (Elbakidze et al. 2013). Also, interactions between plants and the biotic and abiotic properties of soils play an important role in structuring terrestrial communities (Van der Putten et al. 2013, 2016). For instance, gum trees (Eucalyptus globulus) are exotic species in Portugal and, because their litter is more complex, persistent and releases toxic compounds that inhibit decomposition (Zhao et al. 2021a, b), forest dominated by this species will have low support for soil biodiversity. In contrast, deciduous tree species, such as oak trees (Quercus spp.) that are native to Portugal, are likely to better support soil communities, as their litter is easier to decompose and a source of nutrients to the food web (Keith et al. 2009). Due to this strong above-belowground interconnection of communities, soil biodiversity and the services it provides are an essential component of forest ecosystem health (Parker 2010). Thus, knowledge about how far we can impact soils without losing their functions and services is crucial. Habitat destruction, intensified and altered land-use, fragmentation, and increasing population density seem to be key driving forces for ecosystem change and should be minimised to increase forest ecosystem resilience (Leberger et al. 2020; Milad et al. 2011).

Soils regulate the diversity and functioning of aboveground systems, and are particularly important for the decomposition of organic matter and carbon sequestration in forest areas, as well as their benefits to human well-being (Bardgett and van der Putten 2014; Gessner et al. 2010; Guerra et al. 2021a, b; Ruan et al. 2005). It is not possible to guarantee future aboveground biodiversity and food supply if we do not safeguard soils and their biodiversity for future generations (Guerra et al. 2021a, b). Since soil biota are sensitive to changes in their environment, they are often used for early diagnosis of multiple ecosystem function variations (Havlicek 2012). Among them, soil nematodes are widely-used indicators of soil disturbance, quality and functioning because they are ubiquitous and abundant, and have a short response time relative to environmental impacts (Bongers and Bongers 1998; van den Hoogen et al. 2019).

Nematodes are found wherever decomposition occurs and have key roles in regulating carbon fate and controlling soil microbial populations (Bongers and Bongers 1998). Their community composition reflects soil texture, climate, biogeography, organic inputs, and natural and anthropogenic disturbances (Bongers and Bongers 1998; Brussaard 1997; van den Hoogen et al. 2019; Yeates 2003). Nematodes can be classified in functional groups based on their morphology and feeding habits, with soils including complete trophic food webs composed by bacterivores, fungivores, herbivores, omnivores and predators (Yeates et al. 1993). The application of nematode indices that integrate information on the assemblage of nematode functional groups is a useful tool in soil studies. However, a more recent expansion on the indices is the calculation of metabolic footprints that directly indicate the status of soil functions and, indirectly, the provision of ecosystem services (Du Preez et al. 2022; Ferris 2010; Sanchez-Moreno and Ferris 2018). The integrated analysis of nematode communities is a useful tool that can be applied in comparative studies to infer on belowground biodiversity and function in all ecosystems. The conservation status of most soil organisms is seldom addressed, and only recently have the state of soil biodiversity been assessed in response to conservation strategies (Ciobanu et al. 2019; Guerra et al. 2022; Zeiss et al. 2022). Nematode metabolic footprint analysis has been employed to assess the biodiversity and functions of forest soils (Silva et al. 2021; Zhang et al. 2015). To the best of our knowledge, studies that compare soil nematode communities areas according to forest protection status are unavailable.

This study aims to investigate whether and to what extent protected areas and forest typology (dominance of exotic or native trees) affect soil nematode biodiversity and functions, such as decomposition/mineralization of nutrients and pest and disease regulation (Ferris 2010; Sanchez-Moreno and Ferris 2018), in forest ecosystems. We studied soil nematode communities to: (i) assess differences in their biodiversity and functions between protected and non-protected areas between native and exotic forest types; and (ii) infer how protected areas deal with pressures and drivers to soil biodiversity and functions (such as pH, soil moisture and apparent density).

Material and methods

Study area

The North of Portugal is characterised by a remarkable diversity of forest ecosystems, both natural (mainly dominated by oak species) and anthropogenic (planted forests for silviculture). Our study was located in an area of approximately 21,515 km2, with 15.5% occupied by forests that include some of the rarest and most valuable elements of regional and national biological diversity (Tereso et al. 2011). Native forests are represented by common oak (Quercus robur) and Pyrenean oak (Quercus pyrenaica), accompanied by cork oak (Quercus suber), and also pine trees (maritime pine Pinus pinaster and Scots pine Pinus sylvestris). Exotic forests are represented by gum tree plantations (Eucalyptus globulus), and invasive species of the genus Acacia, namely the silver wattle (Acacia dealbata), blackwood (Acacia melanoxylon) and golden wattle (Acacia longifolia) (Carmo et al. 2011; Tereso et al. 2011).

This region is subject to a wide variety of pressures, such as invasive species, land-use change (Carmo et al. 2011; De la Fuente and Beck 2018; Vicente et al. 2013), and forest fires (Alcasena et al. 2021; Carmo et al. 2011; Marcos et al. 2019). Moreover, about 25% of this territory is classified as areas of conservation value, including several areas integrated in the Natura 2000 network, such as Montesinho Natural Park, Douro International Natural Park, Litoral Norte Natural Park, Alvão Natural Park, Tua Valley, Douro Estuary, Serras do Porto Park, and Corno do Bico, and the only National Park in Portugal, Peneda-Gerês National Park (ICNF 2017; Vicente et al. 2013).

Sampling and data collection

In order to have an adequate representation of the local environmental and social-economic conditions in the North of Portugal, we implemented a stratified random sampling scheme considering a 20 × 20 km grid across the entire region. For each grid cell, two sampling sites were randomly selected, one of each forest type [n = 127 samples; subdivided into those dominated by native (n = 70 samples) or exotic (n = 57 samples) trees], the same sites were also classified as inside or outside protected areas (protected = 38 samples, non-protected = 89 samples) (Fig. 1). The samples were taken in the Autumn, from October to December 2021.

Fig. 1
figure 1

Effects of protected areas on soil functional biodiversity in the North of Portugal. a Distribution of the sampling sites (n = 127) within protected areas (n = 38) and environmentally similar non-protected areas (n = 89); b Box plot representing the standardized edaphoclimatic variables used for comparison between sampling sites. *p ≤ 0.05; **p ≤ 0.01; ****p ≤ 0.0001

At each site, a homogeneous soil sample composed of nine soil cores distributed in a 15 m grid pattern was collected from a 30 × 30 m square. Soils were collected down to 10 cm depth using a hori-hori knife after removing superficial plant litter. The soil cores were combined and mixed carefully in a plastic bag, further removing large stones and plant material. A 120 ml plastic container was filled with the sampled soil, transported in a refrigerated box to the laboratory, and kept at 5  ºC for nematode extraction. The remaining disturbed dry soil sample, weighing about 1 kg, was stored at room temperature, approximately 23  °C, for physical and chemical analysis. An additional undisturbed soil cylinder was collected at the centre of the square to assess apparent density. Environmental and edaphoclimatic variables include, land-use classification, mean annual temperature and precipitation seasonality (Karger et al. 2017), distance to urbanisation and distance to rivers [Copernicus Land Monitoring Service 2018, European Environment Agency (EEA)], and soil organic carbon (FAO 2022) (a full list of variables is given in Supplementary material Table SI).

Quantification of the physico-chemical parameters

For each sampling point, a set of physico-chemical parameters, namely apparent density, electric conductivity, magnesium, phosphorus, potassium, soil organic carbon, moisture and pH, was collected. The procedures for collecting undisturbed samples for the determination of apparent density and porosity followed ISO 11272 and 11274 standards. Organic matter/organic carbon was determined by digestion with sodium dichromate and ultraviolet/visible molecular absorption spectrophotometry (Nelson and Sommers 1996; FAO 2020). The pH was determined potentiometrically in a soil:water suspension in the ratio 1:2.5 (m/V). After manually shaking a 1:5 deionised aqueous solution for two hours and then filtering, the electric conductivity value of fresh samples was quantified using an electrical conductivity meter. Phosphorus (P2O5) and Potassium (K2O) contents were determined using the modified Egner−Riehm method (Egner and Riehm 1960). Phosphorus and potassium were extracted from a solution of ammonium lactate and acetic acid measured using colorimetry and flame photometry, respectively. Extractable calcium and magnesium were quantified by atomic absorption spectrophotometry in ammonium acetate (Póvoas and Barral 1992).

Nematode extraction and identification

Nematodes were extracted for 72 h from 100 ml of soil from each sample using the tray method, an adaptation of the Baermann Funnel (Whitehead and Hemming 1965). The obtained ca. 300 ml nematode suspension was concentrated by sieving through a 20 µm-pore sieve. All nematodes were counted in a nematode counting dish, and the first observed 100 identified. When nematode densities exceeded 500 in a given sample, up to 200 were identified. Nematodes were observed using an inverted microscope (Leica DMi1), from 50× to 200×, classified in trophic groups through their mouth and upper digestive tract anatomy (Yeates et al 1993), and to family or up to genus level using simplified keys (Ferris 2011; Goodey 1963; Mai and Mullin 1996; Tarjan et al. 1977). Finally, nematodes were classified in functional groups (Bongers 1990; Yeates et al. 1993).

Statistical analysis

Nematode taxa abundance was pre-analysed using the NINJA (Nematode Indicator Joint Analysis) software (Sieriebriennikov et al. 2014) to calculate the Enrichment Index and the Structure Index (Ferris et al. 2001), and nematode metabolic footprints (Ferris 2010). The latter included the Enrichment and Structure footprints, and that of each trophic group (i.e., bacterivores, fungivores, herbivores, omnivores, and predators), which were plotted on radar charts (Sanchez-Moreno and Ferris 2018). The Enrichment Index reflects the occurrence of opportunistic nematodes (bacterivorous and fungivorous), indicative of food resource enrichment (fungivorous for more recalcitrant sources and bacterivorous for more nutrient-rich sources). The Structure Index is based on the relative abundance of nematodes in higher trophic groups and colonizer-persister levels (Bongers and Bongers 1998), being indicative of soil food web structure and connectance (Ferris et al. 2001). In turn, the metabolic footprints measure the magnitude of ecosystem functions performed by each nematode functional group, with the Enrichment, Bacterivore and Fungivore footprints assessing the decomposition service, the Structure, Omnivore and Predator footprints indicating the natural regulation service and finally, the Herbivore footprint assessing the herbivory disservice (Ferris 2010). The nematode trophic groups abundance and nematode alpha diversity (richness) was calculated based on the number of nematodes per 100 g of soil and the number of nematode families or genera per sampled site, respectively. Variables were transformed using logarithm to normalize the data.

The analyses on the effects of land use and protection on soil nematode abundance, diversity and footprints were performed in R version 4.0.3 (R Core Team 2020). Due to the mismatch on the coupling of samples by grid, differing number of samples for each forest and protection status and heterogeneity of the area, pairs of sites including native and exotic or protected and non-protected forest areas were randomly selected, assigning each native site to an environmentally similar exotic site or, in the case of protection status, each protected site to an environmentally similar non-protected site in 1000 independent runs (total sites n = 127, exotic = 57, native = 70, protected = 38, non-protected = 89). The environmental similarity was quantified as the lowest Mahalanobis distance using the set of variables: soil organic carbon, precipitation seasonality, electric conductivity, soil moisture, magnesium, phosphorous, distance to urbanisation, apparent density, distance to rivers and mean annual temperature (Smith et al. 2021; Zeiss et al. 2022). To exclude highly correlated environmental variables, the Variance Inflation Factor (VIF) was used (Fox and Weisberg 2019). Variables were considered highly correlated if the VIF value was above 5 and eliminated one by one, thus obtaining a shortlist of 24 variables. This was further shortlisted through the analysis of Spearman correlations between variable pairs. The ones that were significantly correlated (p < 0.05) and with a correlation coefficient larger than 0.5 or smaller than − 0.5, were omitted, thus obtaining a final shortlist of 13 environmental and edaphoclimatic variables. Soil nematode biodiversity and functions of all pairs were considered to significantly differ if the mean of p-values of all 1000 runs from the Welch Two Sample t-test was lower than 0.05 (Zeiss et al. 2022).

Similarity in nematode communities between land use types and protection status was assessed by means of Bray–Curtis indexes. We performed a Principal Coordinate Analysis (PCoA) based on the Bray–Curtis dissimilarity index, to check for possible differences in nematode community composition between exotic and native forests, inside or outside protected areas. Furthermore, a permutational multivariate analysis of variance (PERMANOVA) was carried out to see if nematode community composition differed among the study sites; 9999 permutations were carried out using the adonis function in vegan (Oksanen et al. 2011). Normalized Mantel tests, based on 9999 permutations, were performed to test the correlation between Bray–Curtis distances and geographic distance, using the statistical package vegan (Oksanen et al. 2011). Partial Mantel tests (Legendre and Legendre 1998) were used to estimate the correlation between nematode community similarity and geographic distance while controlling for the effect of soil physical chemistry and environmental conditions.

The influence of the environmental and edaphoclimatic variables, selected with the previous criteria, on soil nematode trophic groups abundance, was determined by Canonical Correspondence Analysis (CCA) using the PAST software (Hammer et al. 2001). Significance of the CCA axes was assessed through permutation analysis and the contribution to the ordination of values related to the abundance of nematode trophic groups was assessed by the obtained correlation values.

Results

Effect of protected areas and forest type on nematode trophic abundance and diversity

Soils in protected areas had a greater abundance of fungivorous and omnivorous nematodes than soils in the non-protected areas (Fungivores: Mean protected = 1723.38; Mean unprotected = 1261.25; p = 0.04, 95% CI [0.06; 0.01]; Omnivores: Mean protected = 136.99; Mean unprotected = 70.64; p = 0.03, 95% CI [0.03; 0.01]) (Fig. 2a). This was not observed for the abundance of other nematode trophic groups or nematode alpha diversity (richness), whose values did not differ statistically between the two types of protection status, nor according to forest type. However, fungivorous and omnivorous nematodes were also significantly more abundant in native forest than in exotic forests (Fungivores: Mean native = 1746.66; Mean exotic = 973.22; p = 0.001, 95% CI [0.001; 0.0006]; Omnivores: Mean native = 112.93; Mean exotic = 62.93; p = 0.01, 95% CI [0.02; 0.009]) (Fig. 2b). Nematode community data is presented in supplementary material Table SII , and nematode parameters data, such as trophic groups abundance, nematode diversity and nematode footprints are shown in supplementary material Table SIII.

Fig. 2
figure 2

Box plot representing the abundances of the nematode trophic groups and alpha diversity in a forest areas inside or outside protected areas, and b forest areas containing native or exotic tree species. *Significant difference (p < 0.05). Standardized values: Variables values transformed with logarithm for a normal distribution

The PCoA plot showed that nematode communities were significantly clustered among exotic and native forests (R2 = 0.075, PERMANOVA permutation 9999 times, P < 0.001) and among protected and non-protected forest areas (R2 = 0.022, PERMANOVA permutation 9999 times, P < 0.05) (Fig. 3). Dissimilarity of the nematode communities was positively correlated with geographic distance (r = 0.07, p < 0.05). Similarity decays with geographic distance, however, despite its significance; this correlation is very small, as indicated by the Mantel test. To determine whether this pattern of increasing dissimilarity with distance was caused by changing soil properties we used partial Mantel Test, that showed no significant relationship between nematode community composition and soil properties (r = − 0.02, p = 0.67).

Fig. 3
figure 3

Principal Coordinate Analysis based on Bray–Curtis dissimilarity of soil nematode community structures in a non-protected and b protected forest areas sites, both separately for exotic and native forests

Effect of protected areas and forest type on nematode footprints

Protected soils had a significantly larger Structure Footprint than non-protected areas [Structure Footprint: Transformed mean protected = 2.31; Transformed mean unprotected = 2.11; p = 0.04, 95% CI (0.05; 0.01)] (Fig. 4a). The same was observed for native forests that also had a greater Structure Footprint than exotic forests (Structure Footprint: Transformed mean native = 2.34; Transformed mean exotic = 1.96; p = 0.002, 95% CI [0.003; 0.001]) (Fig. 4b). Irrespective of protection status, forests seemed to have larger Fungivore and Omnivore footprints than Herbivory footprint. The Omnivore footprint showed a significant difference between the two types of protection (Omnivore Footprint: Transformed mean protected = 2.09; Transformed mean unprotected = 1.91; p = 0.04, 95% CI [0.04; 0.007]), with forest soil in protected areas having a larger Omnivore Footprint (Fig. 4a). Nutrient allocation in soil food webs both in native and in exotic forests was larger in the decomposer food web (i.e., larger metabolic footprints of fungivores and omnivores) than in the primary production channel (smaller Herbivory Footprint). In addition, native forests had overall greater nematode metabolic footprints than exotic forests, except for the Bacterivore, Predator and Herbivore footprints, for which significant differences were not detected. The metabolic footprint of fungivores was significantly greater in forest areas dominated by native tree species than in those dominated by exotic tree species (Fungivore Footprint: Transformed mean native = 2.12; Transformed mean exotic = 1.81; p = 0.0009, 95% CI [0.001; 0.0005]) (Fig. 4b). The Omnivore Footprint was significantly greater in native forest areas than in exotic forest areas (Omnivore Footprint: Transformed mean native = 2.10; Transformed mean exotic = 1.80; p = 0.01, 95% CI [0.01; 0.007]) (Fig. 4b).

Fig. 4
figure 4

Nematode metabolic footprints as indicators of the magnitude of the ecosystem functions performed by each nematode functional group for a forest areas inside or outside protected areas, and b forest areas with native or exotic tree species. *Significant difference (p < 0.05) from the Welch Two Sample t-test. Values were transformed using logarithm for a normal distribution

Influence of the environmental and edaphoclimatic variables

Canonical Correspondence Analysis (CCA) showed no relationship between environmental and edaphoclimatic variables with nematode metabolic footprints nor with nematode taxonomic diversity. However, the environmental variables explained 17% of the variance in nematode trophic diversity with 12% and 3% explained on axis 1 and axis 2, respectively (Fig. 5 CCA). Considering the significant axis 1 (p = 0.02) predatory nematodes were clearly separated from the other trophic groups. Both predators and omnivores, that are nematodes closely linked to complex and long food chains, were on the positive side of the axis, thereby positively associated with soil organic carbon and moisture and negatively associated with apparent density. The bacterivores, fungivores and herbivores were on the negative side of the axis, being negatively correlated with soil organic carbon and moisture and positively associated with apparent density. Furthermore, bacterivores were closer to the origin of the coordinate system, being thereby less correlated to the environmental variables. Values of the environmental variables and physico-chemical parameters used in the CCA are shown in supplementary material Table SIV .

Fig. 5
figure 5

Canonical correspondence analysis bi-plot of the influence of the environmental and edaphoclimatic variables on nematode trophic groups. Environmental variables are marked by arrows. Eigenvalues were 0.008704 and 0.002156 for the first and the second axis, respectively (Triplot amplitude: ×3). Forest_Dec deciduous forest cover; SOC soil organic carbon; Forest_Con coniferous forest cover; MAP_Seas mean annual seasonal precipitation; CE electric conductivity; Moist_vol soil moisture; Mg magnesium; P phosphorous; Dist_Urban distance to urbanization; App_Dens apparent density; Dist_River distance to rivers; MAT mean annual temperature

Discussion

When addressing the effect of nature protection in soil nematode biodiversity and functions, the Structure and Omnivore footprints were found to be significantly different. However, we found no evidence for differences between the other nematode metabolic footprints in protected compared to non-protected forest areas in the North of Portugal. The greater Structure and Omnivore footprints in protected areas suggests that these areas maintain a more stable and complex community structure, with a larger investment in natural regulation (Ferris et al. 2001; Gao et al. 2020; Tomar and Ahmad 2009). The nematode alpha diversity did not significantly differ between the protection status, however, since the Structure and Omnivore footprints were significantly greater in protected areas, a tendency towards positive effects of protected areas in the soil system was apparent. Studies of the soil biodiversity and function in protected forest areas are scarce and show limited to no effects of nature conservation on soil biodiversity (Ciobanu et al. 2019; Zeiss et al. 2022). Nematode analyses of Romanian grasslands assigned to Natura 2000, for example, showed no evidence of significant positive effects of protected areas on soil biodiversity and functions (Ciobanu et al. 2019); contrary to what we found here, where protected forests showed some significantly positive effects.

Soil nematode functions were affected by forest type. The Fungivore, Omnivore and Structure footprints were significantly higher in native than in exotic forest areas in the North of Portugal. Native forests, unlike exotic forests, were shown to have a more complex community structure with more linkages in the food web (Ferris et al. 2001; Gao et al. 2020; Tomar and Ahmad 2009). Low abundances of free-living nematodes may indicate lower organic matter decomposition and nutrient turnover rates in gum tree plantations in comparison to natural forests. Exotic forest areas have also been found to have lower nematode abundances in contrast to native forest areas that have higher abundances (Gao et al. 2020). The lower nematode abundances in exotic tree species may be because gum tree litter is more persistent and complex to decompose, and releases toxic compounds that inhibit decomposition, and reduces nutrient availability, supporting fewer nematodes (Zhao et al. 2021a, b).

Studies regarding nematode community dissimilarity have already been carried out for regional scales and in forest ecosystems, where differences were found between forest types (Kitagami et al. 2022), and nematode community similarity decayed significantly with geographic distance (Monroy et al. 2012). Taking into account the PCoA plot, this supports our results where nematode communities were significantly clustered among exotic and native forests and among protected and non-protected forest areas, meaning that nematode communities differed between the two types of forests and protection status. Our results showed that, nematode communities dissimilarity was positively correlated with geographic distance, and that the similarity decays with geographic distance, however, despite its significance; this correlation was very small. This goes in accordance with previous studies where the same was stated (Monroy et al. 2012). According to our partial mantel Test, this correlation is highly dependent on soil properties, as when these variables are taken into account, correlation between community and geographic distance is no longer significant.

Regarding the environmental factors that may influence nematode community structure, vegetation and litter inputs appear to be the most influential factors. Aboveground vegetation can affect the soil nematode community structure, and therefore the diversity and composition of plant species can be an important factor in determining the structure of nematode communities (Kitagami et al. 2017). The quality of litter inputs is a strong determinant of soil community structure. The diversity and composition of the aboveground plant species were not determined in this study, but our results suggest the percentage of deciduous forest cover may have influenced nematode community structure and functions.

The variables that mostly influenced or overall affected the soil functional biodiversity were soil organic carbon, apparent density, soil moisture, and deciduous forest cover. However, they together had a low percentage of variance explained. The influence of these variables was already documented in various studies considering forest areas (Karuri 2021; Keith et al. 2009; Zhang et al. 2012), and these low percentages of variance explanation were also documented by Archidona-Yuste et al. (2020). However, forest areas studies show higher explanation percentages, approximately 35% (Salamon and Wolters 2009; Zhang et al. 2012), but still variables explain only in part the nematode community structure (Salamon and Wolters 2009). An explanation for the lower percentages of variance explanations of environmental variables analysed may be that a large proportion of the unexplained variation is driven by features and mechanisms that were not considered, for example species interactions (Archidona-Yuste et al. 2020), microsite characteristics and land use management (for instance exotic forests tend to be planted, which involves ploughing and terracing soil more frequently than native forest). Thus, effects on microhabitat structure (Bhusal et al. 2014), temperature, soil structure (Salamon and Wolters 2009), elevation (Kergunteuil et al. 2016) or microbial biomass were previously found to be among the main factors determining nematode genera distribution (Zhang et al. 2012). We did not have access to these variables for this analysis, but we expect that they would help to understand the observed patterns in nematode communities. Further studies addressing floristic diversity and elevation, for example, are also needed to complement datasets and possibly increase the proportion of explained variance in this study. This study attempts to interpret macro-ecological patterns at the regional level, and some factors such as plant community structure and land-use intensity can affect nematode community structures on a regional scale (Liu et al. 2019). Other influences such as habitat destruction, fragmentation, and increasing population density may also be key driving forces for ecosystem change and should be minimized to increase forest ecosystem resilience (Leberger et al. 2020; Milad et al. 2011). Promoting connectivity of forest landscapes and habitat heterogeneity is likely to promote species survival, as well as changes in their composition and hence forest management (Milad et al. 2011).

In regard to our first hypothesis, we infer that current conservation measures are having few significant effects on these communities. Furthermore, native forests appear to be more important as a conservation target than exotic ones. Concerning the second objective of this work, we conclude that current conservation measures still need to be improved to assure and increase soil biodiversity and function. It is worth mentioning that protected areas in Portugal include human settlements living as socio-ecological systems, and only very small areas have restricted activities to forest exploitation, therefore a small influence of protected areas on nematode abundance and metabolic footprint was expected. Future research still needs to address species interactions, and to estimate effects at a local scale, e.g., at a case-by-case analysis.

Planting deciduous native trees as mixed forest species has been reported to improve adaptation to site conditions and conservation issues (Felton et al. 2016, 2010; Salamon and Wolters 2009). However, the extent to which ecological, economic, and societal benefits may rise from the conversion of monocultures to mixed forests depends on the tree species being implemented and the bio-geographical region (Felton et al. 2010). It is not clear whether turning the gum tree and wattle monocultures into mixed forests using pine trees would be efficient to improve their environments, as they produce highly recalcitrant leaf litter (Kondratow et al. 2019). However, given our results, we argue that interplanting oak trees in exotic and pine tree dominated forests, at a density and distribution to be determined may possibly increase the ecological and conservation value of these forests. After the tragic forest wildfires in 2017, some efforts are being done to create incentives for a native and more resilient forest cover in Portugal (Santos et al. 2019).

Conclusion

Current work done by nature conservation in the North of Portugal is producing significant effects on the soil system. This highlights the importance of advising nature conservation policies on the potential positive effects on soil biodiversity and functions. Furthermore, native forest areas are better suited to ensure soil biodiversity and functions and should have a greater importance in nature conservation measures. Our study also showed that soil organic carbon, soil apparent density, soil moisture, and deciduous forest cover were the variables that mainly influenced or affected the soil functional biodiversity of nematodes. However, further research is needed to complement, corroborate, or even contradict some conjectures presented in this study, such as the potential effects of species interactions, floristic diversity, and local, microsite-level, and heterogeneity. Here, we concluded that when establishing conservation areas, local soil biodiversity and soil monitoring should be considered. We present evidence for a better belowground forest conservation suggesting that turning current monocultures of exotic species and pine trees into mixed forests may possibly imply significant positive changes to the preservation of soil communities. While such practices have already been successfully implemented in other regions, for the first time, we show their potential consequences for the conservation of soil biodiversity and function.