Changes in Soil Hydro-Physical Properties and SOM Due to Pine Afforestation and Grazing in Andean Environments Cannot Be Generalized

Andean ecosystems provide important ecosystem services including streamflow regulation and carbon sequestration, services that are controlled by the water retention properties of the soils. Even though these soils have been historically altered by pine afforestation and grazing, little research has been dedicated to the assessment of such impacts at local or regional scales. To partially fill this knowledge gap, we present an evaluation of the impacts of pine plantations and grazing on the soil hydro-physical properties and soil organic matter (SOM) of high montane forests and páramo in southern Ecuador, at elevations varying between 2705 and 3766 m a.s.l. In total, seven study sites were selected and each one was parceled into undisturbed and altered plots with pine plantation and grazing. Soil properties were characterized at two depths, 0–10 and 10–25 cm, and differences in soil parameters between undisturbed and disturbed plots were analyzed versus factors such as ecosystem type, sampling depth, soil type, elevation, and past/present land management. The main soil properties affected by land use change are the saturated hydraulic conductivity (Ksat), the water retention capacity (pF 0 to 2.52), and SOM. The impacts of pine afforestation are dependent on sampling depth, ecosystem type, plantation characteristics, and previous land use, while the impacts of grazing are primarily dependent on sampling depth and land use management (grazing intensity and tilling activities). The site-specific nature of the found relations suggests that extension of findings in response to changes in land use in montane Andean ecosystems is risky; therefore, future evaluations of the impact of land use change on soil parameters should take into consideration that responses are or can be site specific.


Introduction
High montane forests and páramo [1], typical Andean ecosystems, provide important services such as water supply regulation and carbon storage [2]. The hydrological services provided by Andean ecosystems have been closely related to the characteristic functions played by their soils at watershed scale [3][4][5]. These characteristics include: high soil organic matter (SOM) content [6], overall low bulk density (BD) (<0.6 g cm −3 ) [7], and high-water retention capacity [7][8][9][10][11][12]. Despite the important role soils play with respect to the provision of ecosystem services in the Andean region, they have been increasingly altered by anthropogenic activities since pre-Hispanic times [13], and only a limited number of evaluations of the effects of anthropogenic land-use changes on soil properties exist.
In the Andean region, afforestation and grazing are the main anthropogenic activities affecting montane ecosystems [12,[14][15][16]. Afforestation, primarily with pine, has been carried out to produce timber and reduce pressure on native ecosystems [16][17][18][19], as well as to capture carbon [19][20][21]. Without focus on the possible impact on soil and water resources, the majority of pine plantations were established on sites of important hydrological and ecological value to the Andean population. Grazing in the montane Andean ecosystems has been carried out since pre-Hispanic times (XV century), and this practice results mainly from extensive to intensive animal farming (cattle, horses, and sheep) [13]. In some cases, grazing is carried out directly in natural cover, considering tussock grass as fodder, or after burning to convert the land, whether or not accompanied with the introduction of more productive grass species [22], into suitable pasture where cattle can obtain more nutritious fodder [18].
Pine plantations and grazing interventions in Andean ecosystems have led to significant changes in the hydro-physical properties and SOM content, in particular of the top soil at depths lower than 30 cm [19,20,23,24]. Modifications in soil properties are causal factors of a reduction in water yield and variations in the hydrological response of montane Andean basins [25,26]. In the case of pine afforestation, the soils tend to dry out because of the high water absorption capacity of the trees, which favors SOM decomposition [19,20] and reduces the soil's water retention capacity [15]. Furthermore, preferential flows created by root growth increase the saturated hydraulic conductivity (K sat ) of soils [27]. Despite existing investigations, there still exists indistinctness and contradictions about the effect of anthropogenic impacts, such as pine plantations and grazing, on the soil properties of Andean ecosystems [28]. Some studies reported an increment in SOM in the Andosols of pine plantations [28], while other studies found no evidence of changes [29]. Additionally, in most studies, important factors capable of influencing changes in soil hydro-physical properties, such as previous land-use practices as evidenced by La Manna et al. [30], have not been considered. Similarly, the effects of livestock grazing, whether or not accompanied by burning, have been interpreted in different ways. Some authors claim that greater solar radiation exposure leads to soil drought, an increase in BD, SOM loss [19,23,31], and a reduction in soil water retention [23,32]; while a reduction of K sat enhances water erosion [28]. Despite previous findings, under grazing activities some authors have shown, probably due to pre-tilling activities, a reduction in the BD and an increment in SOM [28,33]. Such discrepancies could be due to the intensity of grazing. In some cases, extensive grazing is not considered a source of stress on high Andean ecosystems [34].
Considering the role of soil properties in the provision of ecosystem services in Ecuador's Andean environments, and the past and present anthropogenic pressures resulting from afforestation and grazing activities, the main objective of this study was to evaluate the impact of pine afforestation and grazing on the hydro-physical properties and SOM content of the surface horizons of high montane forest and páramo soils. The study was guided by the following research hypotheses: (1) to properly assess anthropogenic impacts it is recommended to collect information of unaltered and altered plots on neighboring comparable sites to avoid that interpretation of results is biased by local site differences and (2) the impact of land-use change on the hydro-physical soil properties is not unique and often masked by other factors such as antecedent land-use, spatial variability, texture, elevation, climate, among other site-specific factors. The evaluation of both these hypotheses will not only improve the research quality and the capability to generalize research findings, but also the effectiveness of future decision making and conservation policies of high-elevation ecosystems in the Andean region. Climate at the study area is influenced by the Pacific coast regime and the Atlantic air currents [6,39]. Average annual solar radiation varies from 11.35 to 14.5 MJ m −2 day −1 , average daily temperature from 12.5 to 6 • C, and average annual precipitation from 900 to 1600 mm at elevations between 2600 and 4100 m a.s.l. [26,40]. The climate regime is bimodal with a rainy season from December to January and a less rainy season from August to September [41]. Rainfall variation is controlled by elevation and topographic parameters [12,26]. The geological material in southern Ecuador corresponds mainly to the Tarqui, Turi, and Saraguro mountain formations, which date back to the Miocene in the Neogene period. These formations are composed of a wide variety of lithological andesite rocks, ash-flow tuffs, pyroclastic flows and volcanic rocks [42].

Selection and Implementation of Study Sites
Seven study sites (Figure 1) expanding between 2705 and 3766 m a.s.l. were selected at three elevation zones: (1) F1 and F2, located below 3250 m a.s.l. corresponding to the high montane forest ecosystem; (2) F3, F4, and Tg5, located between 3250 and 3550 m a.s.l. representing the transition zone between high montane forest and páramo; and (3) Tg6 and Tg7 located above 3550 m a.s.l. being occupied by páramo.
These sites were selected given that they presented intervened zones afforested with pine and converted to grass land nearby undisturbed natural cover zones, all possessing similar topographic and edaphic conditions. We considered undisturbed natural cover in high montane native forest (F) as those sites covered with tree species having a tilted growth, accompanied by a high density of epiphytes and moss, as well as leaf litter accumulation in the ground [36]. In the páramo (Tg), we considered areas with tussock grass >50 cm high and without the presence of anthropogenic intervention indicator-plants, such as Lachemilla orbiculata (Ruiz & Pav.) Rydb. [8], as undisturbed natural cover. Information was collected from pine plantations and grazing areas through field visits and interviews with the landowners [46]. At each pine plantation, the previous land-use, management, plantation age, and the slope of the land, were registered. In the case of grazing, the previous land-use, pre-tilling, and tilling activities, the slope (Sl) of the land, the age of pastures, as well as the animal load (reported as adult bovine units per hectare, ABU Ha −1 ), were registered. A total of 17 plots with a dimension of 24 × 24 m were implemented over the 3 elevation zones, with the exception of a grazing plot in the F4 transition band between high montane forest and páramo due to the complex topographical situation, very much different from the topography of undisturbed natural cover at this elevation band. The location (geographical coordinates and elevation) and slope (Sl) of each plot, and the density or number of trees per plot (SD) were determined. The diameter at breast height (DBH), height (Ht) and canopy diameter (CD) of each tree were also measured using a tape measure and a hypsometer.

Soil Properties Characterization
A first step in the soil characterization was limited to a qualitative description of the soils in the undisturbed natural cover plots following the Food and Agriculture Organization of the United Nations-guide [47]. The qualitative description focused on the surface horizon and included measurement of the thickness, the number of roots per dm 2 , and the structure and texture by touch. A disturbed soil sample (0.5 kg) was collected in the layer 0-10 and 10-25 cm for physical and chemical analysis. In a second prospection and sampling phase the soil hydro-physical properties and SOM content of the 0-10 and 10-25 cm surface layers of the areas with undisturbed natural cover and the areas with impacts of pine afforestation and grazing were determined. Only these two surface layers of the profile were considered since they represent the zone with highest root density and influence of tilling activities. At each surface layer of each undisturbed natural cover and grazing plot, two undisturbed soil samples (100 cm 3 Kopecky rings) were randomly collected together with a disturbed sample of 0.5 kg. In the case of pine plots, the samples were taken at a distance of 75 and 150 cm from the tree trunks. The tree selection was random with a number of 3 trees per plot. At a similar density, K sat was measured by means of the inverted auger-hole method [48] (3 repetitions per site, data shown represent the average of the repetitions). To avoid the effect on the measurement of the water-filled auger below 10 and 25 cm, respectively, the bottom of the upper and lower top layer, a plastic pipe with closed bottom was inserted in the auger holes after saturation of the soil. This enabled measuring subsequently the saturated horizontal hydraulic conductivity of both surface layers. Given the overall low bulk density of the organic top layers the soil can be considered isotropic, assuming that the vertical and horizontal saturated hydraulic conductivity are similar.

Laboratory Analyses
The disturbed samples taken during the qualitative characterization of the soil in the undisturbed natural plots were dried at room temperature (<30 • C) and passed through a 2-mm sieve. The color, pH, and SOM content were determined; the soil color was determined on wet samples using the Munsell color table. The pH was measured with a potentiometer on a 1:2.5 soil:distilled water solution [49], and the SOM content was determined through ignition of the soil at 410 • C for 16 h [50]. The BD and water retention capacity at pressure heads above field capacity were determined on undisturbed samples and are reported as pF values (or the logarithms of the negative pressure heads) corresponding to pF 0 (saturation point; pressure 1 cm H 2 O), pF 0.5 (pressure 3.1 cm H 2 O), pF 1.5 (31 cm H 2 O) and pF 2.52 (field capacity; pressure 330 cm H 2 O). Water content at saturation was considered as a proxy of porosity. Disturbed sieved soil was used for the measurement of the water retention below field capacity at pF 3.4 (2509 cm H 2 O) and pF 4.2 (wilting point; 15300 cm H 2 O), and the SOM content. The equipment used for the measurement of water retention was composed of sandboxes (pFs 0.5-1.5) [51] and two pressure chambers (low pressure for pF 2.52 and high pressure for pFs 3.4 and 4.2; Soilmoisture Equipment Corp., Goleta, CA, USA) [49]. Gravimetric water contents were transformed into volumetric contents. Gravitational water (GW) and water available for plants (AW) were calculated as the difference between the soil water contents at pF 0 and pF 2.52 and at pF 2.52 and pF 4.2, respectively. BD and pF at saturation were considered as proxies for soil compaction.

Statistical Analyses
For the characterization of the pine plantations and vegetation in the grazing plots and the soils in all three land uses, the median and the 25 and 75 percentiles of the distributions were defined given the lack of normality in most of the datasets as evidenced by the results of the Shapiro-Wilk test (p < 0.05). To define if elevation and differences in site as well as sampling depth have an impact on the interpretation of observations, a Spearman correlation analysis (p < 0.05) was performed between elevation, hydro-physical properties, and SOM content measured respectively at 0-10 and 10-25 cm depth under undisturbed natural cover. Subsequently, the soil properties located at the same elevation and depth were compared among the seven study sites. To this end, we applied non-parametric Kruskal-Wallis test [52] (p < 0.05) given the lack of normality in the datasets. If significant differences were identified, the Nemenyi post hoc test [53] (p < 0.05) for multiple pair comparisons was performed. The Mann-Whitney U test [54] (p < 0.05) was applied to evaluate differences in properties between the depth of 0-10 cm versus 10-25 cm.
The soil characteristics of the 0-10 and 10-25 cm soil layers of the land-use-change affected plots were compared to the same characteristics determined in the undisturbed natural cover plots to assess if the soil properties in high Andean ecosystems were affected by pine afforestation and grazing. In the pine afforested plots, first a screening was made between the properties measured at 75 and 150 cm from the trunk. The datasets were grouped if no significant differences were detected. Similarly, comparisons between the measured characteristics were made for both sampling depths (0-10 cm and 10-25 cm). Comparisons were made using the Mann-Whitney U test (p < 0.05). For a better interpretation of the impacts of afforestation, the relationship between the development of pine plantations and the soil properties was analyzed by means of the Spearman correlation analysis and linear regression (p < 0.05). Those analyses were performed on the variables of all trees with the hydro-physical soil properties and SOM content at both depths. The variables that did not present a normal distribution prior to the linear regression analysis were transformed to a normal distribution by Box-Cox method [55]. An analogue approach was followed to assess the impact of grazing. Per study site, soil characteristics at both depths at undisturbed natural cover and grazing plots were compared using the Mann-Whitney U test (p < 0.05).
To better identify and visualize the impacts of pine afforestation and grazing, forest plot graphs were built based on the differences of the medians of the pine plantations (Xtrat) vs. the undisturbed natural cover (Xcn) [56]. For grazing (Xtrat), the same procedure was followed; the differences were transformed into percentages as described by the following equation [56]: (1) Positive percentages represent an increase in the soil properties due to pine afforestation or grazing, while negative percentages represent a decrease, except for the BD. All statistical analyses were carried out by using the R program, version 3.3.2 [57].

General Description of the Experimental Sites
In the following description of the pine afforested and grazing sites, there is always reference made to the previous land-use and/or native vegetation, as to highlight the changes introduced by respectively the introduction of pine plantation and the occupation of the land by cattle. Table 1 presents the general characteristics of the pine plantations in the different study sites. The F1Pi, F2Pi and Tg6Pi plantations were established on sites where previous land-use comprised native vegetation cover. The slope on those sites varied from 43 to 20, and 22%, respectively. The F3Pi, F4Pi, and Tg7Pi plantations were established on sites where extensive-grazing had taken place for more than 5 years with an animal load <1 ABU Ha −1 . These plantations had a slope of 12%, 16%, 34% and 20%, respectively. The majority of the plantations were 16 to 19 years old, except in the F1Pi and Tg5Pi sites where the age of the trees varied between 20 and 29 years. The number of trees per plot (SD refers to the number of trees per 576 m 2 ) varied from 30 to 49 trees. The SD on the F1Pi, Tg6Pi, and Tg7Pi sites was smaller and varied between 30 and 34 trees. F2Pi, F3Pi, F4Pi and Tg5Pi had a greater SD (40-49 trees) (Table 1). Tree development was greater at F1Pi, F2Pi and Tg7Pi with a diameter at breast height (DBH) >18 cm, a tree height (Ht) >9 m, and a canopy diameter (CD) ≥5 m. The pine trees had slower growth in the Tg5Pi, Tg6Pi and Tg7Pi plantations, with an average DBH <18 cm, a Ht <8.5 m and a CD <5 m. In general, the management of the plantations was deficient because of the high associated costs; the F1Pi, F2Pi, Tg5Pi and Tg7Pi plantations did not receive any kind of management, while the only intervention in the F3Pi, F4Pi and Tg6Pi sites was the sporadic pruning of the trees. Table 2 provides a general description of the grazing activity on each site controlled by cattle farming. Grazing takes place primarily on established sites with introduced grasses (G) (F1G, F2G, F3G, Tg6G and Tg7G) except on the Tg5G * site where the cattle grazed the native cover of tussock grass (G *). In all cases, the previous land use corresponded to native vegetation. The type and intensity of pre-tillage and management tasks were different at each site ( Table 2). The age that land was converted to pasture was greater than 10 years in F1G and F2G, unlike F3G, Tg6G and Tg7G where it was shorter, between 3 and 7 years. In Tg5G *, the grassland renewal time was <3 years because of burning. Extensive grazing took place at the Tg5G *, Tg6G and Tg7G sites, with an animal load between 0.5 and <0.2 ABU Ha −1 , while in F1G, F2G and F3G grazing was intensive with an animal load between 1 and 2 ABU Ha −1 .

Soil Properties in Undisturbed Natural Cover Sites
The soils in most of the study sites (F2-Tg7) are identified as Andosols (>80%), except in F1 where they were identified as Cambisols. The deepest soils corresponded to the high montane forest sites (0.73->2.28 m), while the shallower soils are in the páramo (0.46-1.67 m). The thickness of the surface horizon in the Andosols (Ah) varies from 34 to 106 cm, whereas in the páramo the surface horizon was less thick than 55 cm (Table 3). The root density in this horizon was lower in high montane forest sites (10-64 roots dm −2 ) than in the páramo (30-200 roots per dm −2 ). In both ecosystems, the Ah horizon is black (10YR 1.7/1 to 7.5YR 1.7/1), with a high SOM content (11.38%-42.49%), pH values ranging slightly from acidic to very acidic (6.32 to 4.89), and a structure between granular to block ( Table 3). The texture ranged from loam to clay loam in the high montane forest and from loam to loamy silt in the páramo. The surface horizon in the Cambisols (A) is characterized by a thickness ranging between 36 and 50 cm, low root density (11-33 roots dm −2 ), a brownish color (7.5YR 3/2-10YR 2/2), relatively low SOM content (7.09%-14.75%), and an acidic pH (4.58-5.64). Its texture is intermediate to fine, and the structure ranges from block-like to granular (Table 3). The hydro-physical properties and SOM content of the 0-10 cm and 10-25 cm soil layer, of the 20 study sites are listed in the Appendix A (Table A1 depicts the data for the undisturbed natural  cover areas, Table A2 illustrates the same data for the pine afforested sites, and Table A3 for the pasture sites). Those tables show the median and 25 and 75 percentile value for each of the soil layers and 11 soil parameters evaluated. As revealed in Table A1, depicting for each depth the parameter values for the 7 sites in undisturbed natural cover area, the Kruskal-Wallis test (p < 0.05) showed significant differences between pairs of study sites. Also, the Nemenyi post-hoc test showed significant differences in the characteristics of the 0-10 cm surface layer. For example, K sat values in the F1 and F3 sites in the high montane forest were higher than (12.91 and 17.30 cm h −1 ) and differed significantly from the K sat in the Tg6 and Tg7 sites situated in the páramo (1.38 and 1.92 cm h −1 ). The water retention capacity between pF 0 and pF 2.52 was lower for F1, F2 and F3 (<0.73 cm 3 cm −3 at pF 0 and <0.52 cm 3 cm −3 at pF 2.52), differing significantly from the respective values for Tg6 (0.85 cm 3 cm −3 at pF 0 and 0.66 cm 3 cm −3 at pF 2.52). The highest contents of GW were found for F2 and F3 (0.21 and 0.24 cm 3 cm −3 , respectively), whereas the lowest was identified at Tg5 (0.11 cm 3 cm −3 ). On the other hand, the AW content for F1 and F2 were the lowest (0.08 and 0.11 cm 3 cm −3 , respectively), whereas the highest was for the Tg7 site (0.25 cm 3 cm −3 ). In addition, for the BD, water retention capacity at pF 3.4 and pF 4.2, and SOM content, significant differences were recorded regardless of the type of ecosystem. For example, the water retention capacity at pF 3.4 and pF 4.2 in Tg6 was significantly higher than in Tg7 even though both are located in the páramo.
The differences between study sites are smaller at a depth of 10-25 cm. The high montane forest sites presented higher K sat values than the soil in the páramo plots. Important differences were observed in F1 (2.98 cm h −1 ) against Tg5, Tg6 and Tg7 which are characterized by lower K sat values (<0.46 cm h −1 ). The highest volume of GW was found in F3 (0.24 cm 3 cm −3 ), much different from the low GW volume in Tg5 (0.11 cm 3 cm −3 ). The AW values of the F1 (0.06 cm 3 cm −3 ) and F3 (0.06 cm 3 cm −3 ) sites were significantly lower with respect to the AW values in the Tg5 and Tg6 sites (0.20 and 0.21 cm 3 cm −3 , respectively). Similarly, the water retention capacity in the tension range of pF 0 to pF 3.4, and the SOM content of the sites showed significant differences, regardless of the ecosystem type (Appendix A: Table A1). On the other hand, the water retention capacity at wilting point (pF 4.2) did not show significant differences among the different study sites. Furthermore, at each study site the top layer (0-10 cm) presented higher values for K sat , water retention (pF 0 to 2.52), GW, AW, and SOM as compared to the values of the same properties measured in the second layer (10-25 cm), except for the BD, which increased. According to the Mann-Whitney U test (p < 0.05), K sat differed significantly between the two depths at each study site (Appendix A: Table A1). Likewise, BD increased significantly in the second soil layer at the F1, F3, Tg6 and Tg7 sites. The water retention capacity at pF 0 was greater in the 0-10 cm soil layer and tended to decrease significantly in the layer below (10-25 cm below surface) on the F1, F3 and Tg6 sites ( Table 3). Correlation of the median parameter values and the elevation of the sites, using the Spearman test, showed that most of the parameters and SOM content were strongly correlated with the elevation ( Table 4). The strongest correlations were found for the surface layer (0-10 cm), reflecting a negative correlation of the elevation with K sat (ρ = −0.74, p < 0.05) and a positive correlation with the water retention capacity in the tension range of pF 0 to pF 2.52 (ρ = 0.66-0.71, p < 0.05), as well as with the AW (ρ = 0.73, p < 0.05). The hydro-physical properties and SOM content in the 10-25 cm layer showed weaker correlations with elevation, with the exception of K sat (ρ = −0.66, p < 0.05). Table 4. Spearman correlation coefficients (ρ) between elevation and the hydro-physical properties and SOM content of the soil layer at 0-10 cm and 10-25 cm in the undisturbed natural cover areas.

Changes in Hydro-Physical and SOM Content under Pine Afforestation
The properties found at each of the pine afforested and grazing sites in the 0-10 and 10-25 cm soil layers are presented in the Appendix A, in Table A2 (pine plantation) and Table A3 (grazing). According to the Mann-Whitney U test (p < 0.05), the properties in the pine plantations in the 0-10 cm layer were significantly different from the characteristics measured in the 10-25 cm soil layer. Therefore, the impacts of pine afforestation were evaluated considering both depths. On the other hand, there were no significant differences (p > 0.05) in the properties between sampling distances (75 and 150 cm) from the trees in each of the pine plantations (Table A4), suggesting that sampling could have been carried out at either of these distances. Therefore, the properties' values measured at a distance of 75 and 150 cm from the tree trunk were grouped by depth.
Pine afforestation produced significant changes in the hydro-physical properties and SOM content in both layers as depicted in the forest plots of Figure 2 and Appendix A Table A2. K sat of the upper layer decreased significantly in the F1Pi and F3Pi plantations (61.72% and 79.06%, respectively), which are located in the high montane forest area. In the páramo plantations, the K sat value increased significantly in two of the three plantations (Tg5Pi and Tg6Pi) by 25.23 and 59.22%, respectively (see Figure 2a). The BD value did not show significant changes in most of the sites except for F3Pi which increased significantly by 18.42% (Figure 2b). As further depicted in Figure 2c, the water retention capacity at saturation (pF 0) in the upper layer decreased significantly in the F1Pi, F3Pi, Tg6Pi and Tg7Pi sites, with a reduction between 2.53% and 11.10%, while in F2Pi the water content at saturation increased significantly by 6.30% (Figure 2c). The changing tendency in water retention capacity in the pF range from 0.5 to 2.52 (see Figure 2d-f) was similar to the change in the water retention capacity at pF 0 in most of the plantations, supporting the conclusion that pine afforestation led to significant changes in the soil's hydro-physical properties. The water retention capacity at pF 3.4 and pF 4.2 (Figure 2g,h) decreased meaningfully between 12.21% and 27.52% in the F1Pi, F4Pi, and Tg6Pi sites, but differently in the Tg7Pi site where the increase was significant, varying between 15.54% and 44.31%. The GW in F4Pi, Tg5Pi and Tg6Pi increased significantly (13.15%-39.37%) (Figure 2i), while AW decreased considerably in F3Pi and Tg7Pi, with a reduction varying between 42.50% and 57.64%, (Figure 2j). The SOM content in F1Pi, F2Pi and F3Pi decreased considerably, by 28.97%, 29.00% and 47.90% ( Figure 2k); while in Tg7Pi the SOM content was significantly reduced by 77.62%.
The K sat increased significantly in the 10-25 cm layer in the F1Pi, F2Pi, Tg5Pi and Tg6Pi sites by 84.64%, 30.28%, 305.24% and 140%, respectively; while in F3Pi the K sat was significantly reduced by 52.98% (Figure 2a). There was a significant increase of BD in F3Pi (9.81%), in contrast to Tg7Pi where there was a significant reduction of 2.47% (Figure 2b). Regarding the water retention capacity at saturation, F1Pi declined significantly by 10.62%, contrasting with F2Pi and Tg7Pi where there was a significant increment by 9.60% and 6.32%, respectively (Figure 2c). Similar results were obtained for water retention capacity in pF 0.5 to pF 1.5 (Figure 2d,e). Water retention at field capacity decreased significantly in F1Pi, F3Pi and Tg6Pi by 18.92%, 5.68% and 2.2%, respectively, whereas it increased significantly at F2Pi by 9.06% (Figure 2f). The retention capacities at pFs 3.4 and 4.2 decreased significantly between 2.62% and 30.35% at sites F1Pi and F4Pi (Figure 2g,h); whereas they increased significantly between 3.90 and 25.56% at site Tg7Pi. The GW did not show significant changes in any of the plantations. AW in the F2Pi and F4Pi plantations presented a significant increase of 51.55% and 79%. The SOM content in F1Pi, F2Pi and F3Pi decreased considerably by 55.50%, 9.32% and 23.50%, respectively, differently from Tg7Pi where there was a significant increase by 77.62% (Figure 2k).
Correlation analysis revealed that K sat in the 0-10 cm and 10-25 cm layers presents strong correlations with tree development variables DBH, Ht and CD, and there is clearly an increase of K sat with increasing tree height (data not shown). Likewise, the BD in both layers was also positively correlated with the tree development variables (ρ ≥ 0.31, p < 0.05). That is to say, with an increase in tree growth, there is a significant increment in soil BD. On the other hand, water retention capacity and SOM contents are negatively correlated with tree development variables (ρ ≥ 0.31, p < 0.05). That is, in both soil layers (0-10 and 10-25 cm), water retention capacity and SOM content decrease significantly with tree development.

Changes in Hydro-Physical and SOM Content under Grazing
The Mann-Whitney U test (p < 0.05) revealed significant changes in the soil's hydro-physical properties under grazing, which was reflected mainly in the K sat and water retention capacity (Appendix A: Table A3). There were a higher number of sites with significant changes in the studied soil hydro-physical properties in the 0-10 cm surface layer, than in the 10-25 cm soil layer (Figure 3). In the upper layer, K sat decreased significantly in the F1G and Tg5G * sites by 68.87% and 75.75%, respectively, while there was a significant increase of 55.19% in Tg6G (Figure 3a). With respect to the BD in F2G and F3G, there was a significant increase of 64.29% and 28.25% (Figure 3b). The water retention capacity at saturation in F3G and Tg6G decreased significantly by 10.10% and 5.18%, respectively, while it increased significantly by 8.97% in Tg5G * (Figure 3c). The changes in water retention capacity at pF 0.5 and pF 1.5 were similar to those at pF 0 (Figure 3d,e). Water retention capacity at pF 2.52 decreased significantly in the F3G, Tg6G and Tg7G sites by 13.11%, 8.61% and 9.77%, respectively. On the other hand, there was a significant increment of 10.46% at F2G (Figure 3f). The water retention capacity at pF 3.4 increased significantly at the sites F2G, Tg5G * and Tg7G by 17.38%, 6.15% and 20.46%, respectively. Conversely, the water retention at field capacity decreased by 7.36% in F3G. The changes in water retention at wilting point (pF 4.2) were similar to those at pF 3.4 (Figure 3g,h) at the sites with significant changes, except for F3G. Changes in GW (Figure 3i) were noticed in the Tg6G and Tg7G sites, with an increase of 10.59% and 36.12%, respectively. The AW (Figure 3j) and the SOM content ( Figure 3k) did not show significant changes in the majority of the study sites, except the for SOM content in site F3G where there was a significant reduction of 30.47%.
In the soil layer 10-25 cm, K sat decreased significantly by 71.01% in the F1G site, by 1120.09% in F3G, by 121.28% in Tg6G, and by 61.35% in the Tg7G site (Figure 3a). The BD did not show significant changes in the grazing sites, except for Tg7G where it decreased significantly by 18.63% (Figure 3b). The water retention capacity at pF 0 increased significantly by 10.71%, 4.67% and 9.57% in the F2G, Tg5G * and Tg7G sites (Figure 3c). At these sites, the trends followed by the changes in water retention capacity at pF 0.5 and pF 1.5 were similar to those at pF 0 (Figure 3d,e). At pF 2.52, the water retention capacity at F2G and Tg5G * increased significantly by 19.17% and 5.92%, respectively (Figure 3f). The water retention capacity at pF 3.4 increased significantly by 12.92% at site F2G, unlike F3G where a significant reduction of 13.78% was found (Figure 3g). Most sites did not display significant changes in the water retention capacity at wilting point, except for F3G, where there was a significant reduction of 12.06% (Figure 3h). Finally, GW, AW, and the SOM content did not show significant changes at all study sites (Figure 3i-k).

427
Some studies in sites close to the study presented herein evaluated the impacts of anthropogenic activities by grouping soils by their conditions of natural undisturbed land cover, without

Discussion
Some studies in sites close to the study presented herein evaluated the impacts of anthropogenic activities by grouping soils by their conditions of natural undisturbed land cover, without considering their elevation and sampling depth [19,28,29]. In our study, results suggest that in order to evaluate the impacts of any anthropogenic activity, it is essential to examine the spatial variability of soils under natural conditions or at least to avoid evaluating the impacts based on sampling sites located in different places. In other words, the impacts of anthropogenic activities on soils should be done in adjacent sites, which would guarantee similar geomorphological and climatic conditions [58]. Even the spatial variability of soils is often the cause that the hydrological performance of Andean ecosystems can be extremely heterogeneous [26,59]. That is, the high variability of soils under natural conditions can even cause variability in the functioning of basins.
An aspect often discussed in monitoring the impact of afforestation on soils is at what distance from the trunk and what depth should samples preferably be collected. A previous study in Andean highlands suggested that changes in soil properties in pine plantations differ depending on the sampling distance from the trunk of the tree [15]. However, in all the plantations within our study area, the hydro-physical properties did not differ between the tested sampling distances of 75 and 150 cm away from the trunk. These results are consistent with what Wilcox, Breshears & Turin [60] and Ruiz et al. [61] reported. Wilcox et al. [60] stated that the K sat of the soil under the canopy and between the canopy of the pine trees did not show significant differences. Likewise, Ruiz et al. [61] reported that the water retention capacity and AW in the soil did not differ between 50 cm from the base of the trunk vs. under the crown of the tree. On the other hand, according to our results, the hydro-physical properties and the SOM content showed significant differences between the two sampling depths (0-10 cm and 10-25 cm) in each of the plantations, a finding in line with the results of Ghimire et al. [62]. Our results suggest that pine afforestation effectively impacts in a different way the soil properties as a function of depth, but not as a function of the distance samples are collected from the tree trunk.
The comparisons of the properties of each plantation with their respective adjacent natural cover revealed that the intervention of the soils through, for example, pine afforestation in the Andean ecosystems mainly affects K sat , water retention capacity between saturation (pF 0) and field capacity (pF 2.52), and SOM content, and this at both sampling depths but at a different intensity. Due to these changes, pine plantations could directly alter the ecosystem services such as water regulation and storage [25,26] and the carbon sequestration by soils [63]. However, findings could not be generalized, showing dependencies mainly in the sampling depth, ecosystem type, characteristics of the plantations, and previous land-use. This implies a complexity in assessing the impacts of plantations and limits the generalization capability of changes in soil properties caused by pine plantations in high Andean ecosystems.
According to Alarcón et al. [64], changes due to grazing were not statistically evident in the K sat measured in the 0-10 cm soil layer in Andosols. However, our results partially contradict this because 2 (F1G and Tg5G *) of the 6 study sites showed a significant reduction in K sat up to 70%, while in the Tg6G site a significant increase of 52.67% was registered. The reduction in K sat at the F1G site could be due to the loss of stability of the soil aggregates and grazing density, while the reduction in K sat at Tg5G * is likely the consequence of the frequent burning of tussock grass resulting in a drying and crusting of the soil surface [65]. On the other hand, although the increase in K sat was not significant at Tg7G, this percentage increase was very similar to that at the Tg6G site (Figure 3a). The increase at both these sites is likely related to preferential flows between the clods formed by soil tillage [66] during the preparation and sowing of pasture. However, this situation could apparently change over time due to the structural deterioration of the soil.
Despite the lack of significant evidence, it was observed that in 6 of the 7 pine plantations, BD tended to increase in the 0-10 cm soil layer, which certainly cannot be the consequence of the use of heavy machinery for maintenance since the associated high costs and the difficult topographic conditions prevent the use of machines [67]. Rather the drying of the soil by evapotranspiration [68] and the weight of the trees [69] are responsible for the increase in BD, which is confirmed in our study by a greater increase in BD in the sites where the pine plantation is characterized by high SD and CD (e.g., F2Pi and F3Pi in high montane forest and Tg5Pi in páramo). The BD at 10-25 cm depth was not affected in most of the plantations; however, most plantations showed a decreasing trend ranging from 2% to 18%. This tendency is attributable to the increase of porosity generated by the pine subsurface root system. On the other hand, two plantations (F3Pi and Tg5Pi) showed an increase which is attributable to the compaction by the pressure exerted by the biomass of the plantation on the soil, due to its high SD, which translates into a greater number of trees, increasing the transpiration, and resulting in soil contraction.
The compaction of the soil by grazing is normally directly reflected in an increase of the BD [70,71]. However, this effect was only observed in two study sites (F2G and F3G) and the direct consequence of the greater grazing intensity at these sites (ABU Ha −1 of 2 and 1, respectively), parallel to a decomposition of SOM. Our findings are consistent with those of Donkor et al. [72] where an increase in compaction was directly related to a greater grazing intensity. In páramo, the trampling effect of cattle was not reflected in the measurements, which was consistent with the results of Alarcón et al. [64] and Podwojewski et al. [23]. Their results and our study suggest that Andosols (due to their high SOM content) have a greater resilience to compaction [73] given the overall low-grazing intensity (ABU Ha −1 of <0.5). Because of the lack of significant evidence, our results do not permit to conclude that grazing significantly compacted the soil at the 10-25 cm depth. This may be due to the high SOM content of the Andosols. In the case of Tg7G, where BD decreased significantly by 18.63%, this decrease is likely caused by the incorporation of tussock grass biomass during tillage.
Pine plantations alter the water retention capacity of soils according to Farley et al. [20]. This is confirmed by the results of this study, with the difference that the alteration is not only dependent on associated changes in the SOM content, but also on elevation, ecosystem type, the development level of the plantation and land-use, as reported in other ecosystems [56,63,74]. The highest reduction percentages in our study were attributed to a greater development of the plantations (DBH >18 cm, Ht >8 m and CD >5 m), together with a decrease in SOM and increase in BD. In the three plantations that were established on soils whose previous use was grazing (F4Pi, Tg5Pi and Tg7Pi), the changes in water retention capacity were very variable, observing significant increases or decreases between pF 0 and pF 2.52, notwithstanding the SOM content tended to increase. This variation in the changes could be due to a mixture of changes between the previous use of the soil and the growth of the plantation [61,74]. This apparent overlap of impacts [62] hinders generalization of changes. On the other hand, our data show an increase in water retention is not necessarily always related to an increase in SOM content, but can also be the consequence of an increase in clay content due to the weathering of the soil [75,76]. In general, pine afforestation goes hand in hand with an increase in GW and a decrease in AW, suggesting that soils under pine plantations rapidly lose moisture after a rainfall event, which further enhances soil drying and decomposition of SOM [77]. On the other hand, according to Buytaert et al. [22], land use change could increase the AW of soils of volcanic origin by 30%. On the contrary, Hofstede et al. [19] reported a decrease in AW as a consequence of pine afforestation. Nevertheless, in our study, only two pine plantations showed a significant decrease (F3Pi and Tg7Pi).
Under grazing, it has been shown that the loss of water retention is mainly a result of a reduction in SOM content. Increases in water retention capacity were explicitly related to increases in SOM content due to burning and/or the incorporation of tussock grass during soil preparation. This finding is in line with Daza et al. [32] who reported that in their study the decrease in soil water retention was due to a loss of SOM content. However, in the F2G site the increase of the water retention capacity was the result of the easy weathering of volcanic glass, leading to the formation of montmorillonite clay [76]. The latter indicates that in order to correctly assess changes in, for example, the hydro-physical properties in soils, it is essential to evaluate the full spectrum of soil properties such as SOM, soil texture, type of management, among others.
Notwithstanding that several studies associate the reduction of SOM content with pine afforestation [19,20,58], our findings indicate that pine plantations could help the recovery of SOM in the 0-10 cm soil layer of former grazing sites, though the effect seems also to depend on the elevation of the site and the SD of the plantation. Although the SOM content did not show significant changes under grazing, a reduction in most of the study sites was detected. This slightly decreasing trend would imply that, over time, the soils under grazing in the Andean region could lose a considerable amount of SOM, and thus, reduce their capacity to retain water.

Conclusions
Andosols with their black Andic horizon are the predominant soils in Andean montane ecosystems. Our research clearly revealed that this horizon, with high water retention capacity and SOM content, is not that uniform under natural unaltered conditions. Differences in hydro-physical properties, such as K sat and the water retention capacity, are related to the type of ecosystem and elevation of the terrain. Differences in properties not only occur between sites, but also within sites at different depths of sampling. As a result, higher values of K sat , water retention capacity and SOM content were recorded in the 0-10 cm surface layer. The natural spatial variability in environmental conditions and the accompanied heterogeneity in soil properties requires that for the correct assessment of the impact of land use change, data are collected on neighboring comparable, unaltered and altered sites. Doing so will facilitate and help guarantee that a correct assessment of the causal factors that positively or negatively affect the soil hydro-physical properties by land use change, is drawn. Furthermore, the multitude of observations and their analyses clearly revealed that the impact of land-use change on the hydro-physical soil properties is not unique and often masked by other factors such as the antecedent land-use, spatial variability, pre-tilling and tilling activities, soil texture, elevation, climate, among other site-specific factors. Due to these differences, it is rather difficult to evaluate the impacts of pine plantations and grazing on properties at regional scale, and therefore any evaluation of the impacts of anthropogenic activities must be carried out in adjacent sites, which would guarantee similar geomorphological and climatic conditions. This conclusion clearly points out that generalization of findings related to the impacts of land use change on soil properties is not free of risks. Similarly, it also hinders the comparison of findings with published results.
The study further revealed that pine afforestation affects either in a positive or negative way the K sat , the water retention capacity in the range pF 0 to 2.52 and the SOM content of the soil surface layer. The change and the order of magnitude of the change varies with sampling depth. Similarly, grazing causes positive and negative changes in K sat and in the water retention capacity, and as in pine plantations, the recorded changes vary with sampling depth. Other controlling factors that define the impact of grazing are evidently pre-tilling and tilling activities in combination with cattle density. Soil spatial heterogeneity and diversity in local factors complicates the interpretation and extrapolation of observed phenomena to Andean montane ecosystems at a regional scale. Correct assessment of land use change impacts is not only of crucial importance for extrapolating findings, it also serves as a basis for the accurate estimation of socio-economic and ecological impacts that anthropogenic-induced changes might have on the water regulation and water storage functionalities of the Andic soils.

Conflicts of Interest:
The authors declare no conflict of interest. Table A1. Median and the 25th and 75th percentile of the hydro-physical properties and SOM content of the soil in the undisturbed natural cover areas (Azuay province, southern Ecuador).       (1) different capital letters represent the significant differences (p < 0.05) between the two depths in each study site (i.e., 0-10 cm versus 10-25 cm). Letter denotation order was used from highest to the lowest value. ↓ or ↑ stand for a significant decrease or increment (p < 0.05 value) when comparing pine plantations vs. undisturbed natural cover.   The names of the variables are defined in Table 4. The values in bold indicate that the properties are significantly different between both depths (p < 0.05).