Shade-Grown Coffee in Colombia Beneﬁts Soil Hydraulic Conductivity

: Secondary tropical forests and coffee agroforestry systems contain fewer trees than native forests but can positively impact soil hydrological functions, such as water inﬁltration compared to the pasture land that they replace. However, for both land uses it remains in how far the soil hydraulic characteristics are comparable to that of native forest. Therefore, we investigated the saturated hydraulic conductivity ( K s ) and some hydrophysical soil attributes in four land-use types: (i) a shade-grown coffee; (ii) a natural regenerated forest 15 years ago; (iii) a pasture; and (iv) a reference forest, in the municipality of La Jagua de Ibirico, C é sar department, Colombia. We determined historical land use and conducted soil sampling, using the Beerkan method to determine the K s values. We also measured canopy cover, vegetation height, diameter at breast height and total number of trees in the forest covers. Our results indicate that K s values were similar for the coffee and the reference forest, reﬂecting the positive effect of trees on soil hydrological functioning in agroforestry systems. Our results suggest that 15 years of forest regeneration after land abandonment in Sub-Andean Forest, can improve the soil hydraulic attributes. Additionally, soil water repellency was observed for the reference forest soil.


Introduction
Agriculture is important for the local livelihoods, the economy, and feeding the world, but it is also a leading cause of environmental damages, including deforestation, biodiversity loss, the spread of invasive species, water pollution and aquifer depletion, soil erosion and land degradation [1,2]. Especially in tropical regions, the use of agroforestry practices has increased over the last decades in order to improve the sustainability of croplands. Agroforestry systems are based on the principle that trees provide benefits to the cropping system [3,4]. In the American tropics, coffee is considered one of the most important economic crops, and is usually produced in agroforestry systems [5]. This is also the case for Colombia, where coffee is the main commercial agricultural crop, generating 12.4% of the total agricultural revenue, moreover, it is estimated that more than half of Colombian coffee is shade-grown [6].
New forests have been emerging in human-impacted landscapes as a consequence of forest restoration actions or natural regeneration [7,8]. These new forests and the shadegrown coffee agroecosystems are important to conserve a large part of the biological diversity of the tropical forests, and also to provide a variety of other ecosystems services (i.e., carbon sequestration, pollination, and production of arboreal biomass) [9,10]. They The study area is located in the municipality of La Jagua de Ibirico, César department, Colombia ( Figure 1). The climate in this region is semi-humid according to the Thornthwaite classification, the elevation is approximately 1300 m, mean annual temperature is 24 • C and mean annual rainfall is 2185 mm [20]. The wettest period is from September to November, while the driest is from December to March [21]. The study area is located in the sub-basin of the Tucuy River, which currently faces several environmental threats including loss of vegetation coverage as a result of mining, agriculture and livestock grazing [22]. The native Sub-Andean forest only remains in remnants of different sizes and ages [20]. The area belongs to a morpho-structural unit known as "Guatapurí", which is characterized by a very pronounced topography [23]. The main soil type in the study area are Entisols (i.e., Typic Udorthents) according to USDA Soil Taxonomy [24]. The soil is shallow and rocky, with a sandy loam texture in the superficial horizon.
Within the area, we selected four land-use types to measure soil physical and hydraulic properties, and vegetation attributes (Figure 2). In each land-use, we established two circular plots, with an area of 79 m 2 each, resulting in eight plots in total. Our study plots were located in a similar landscape position with similar slope (i.e., greater than 20%) and had the same soil type. The size and number of plots were chosen according to previous investigations [14,15].  The CO site (9°35′41.27″ N, 73°6′37.16″ W) is a 35-year-old polyculture with a low use of agrochemicals. The coffee plantation was established with a plant density of 5000 trees per hectare. Our CO site is associated with trees species (i.e., Handroanthus chrysanthus, Inga sp., Persea americana and Tabebuia rosea). The measurements at this site represents the effect of a shade-grown coffee plantation.
The NR site (9°35′52.03″ N, 73°6′54.21″ W) was cleared first in 1960 by slash and burn. After 40 years of cattle grazing the land was abandoned in the early 2000's as a result of the Colombian armed conflict. The area underwent natural reforestation and succession  The CO site (9°35′41.27″ N, 73°6′37.16″ W) is a 35-year-old polyculture with a low use of agrochemicals. The coffee plantation was established with a plant density of 5000 trees per hectare. Our CO site is associated with trees species (i.e., Handroanthus chrysanthus, Inga sp., Persea americana and Tabebuia rosea). The measurements at this site represents the effect of a shade-grown coffee plantation.
The NR site (9°35′52.03″ N, 73°6′54.21″ W) was cleared first in 1960 by slash and burn. After 40 years of cattle grazing the land was abandoned in the early 2000's as a result of the Colombian armed conflict. The area underwent natural reforestation and succession The PA site (9 • 35 46.85 N, 73 • 6 58.63 W), which is the main land-use type in the landscape, was initially a forest that was cleared in the middle of the 20th century. After clearing, the site was heavily grazed for more than 40 years. The area had to be abandoned in the early 2000 s and the vegetation naturally emerged for five years. In 2007 the vegetation was cleared and burned again. Information from landholders revealed that the pasture site has been frequently burned over the last 15 years, but the stocking rate has been low with less than one head of cattle per hectare. At the moment of our field campaign, the PA site was covered by grass species Urochloa brizantha and bracken (Pteridium sp.). Our measurements at this site represent the effect of 15 years of intermittent burning and grazing.
The RF site (9 • 35 48.83 N, 73 • 7 2.22 W) is a secondary old-growth forest that is more than 50-years-old. The RF site was probably burned, similar to the other study sites and there has been selective cutting of trees. Currently, the RF is affected only by the occasional crossing of cattle. We used the RF site as a control area to assess reference values for soil physical and hydraulic attributes.
A graphical summary of the land-use history for the four land-use types described previously is provided in Figure 3. The PA site (9°35′46.85″ N, 73°6′58.63″ W), which is the main land-use type in the landscape, was initially a forest that was cleared in the middle of the 20th century. After clearing, the site was heavily grazed for more than 40 years. The area had to be abandoned in the early 2000′s and the vegetation naturally emerged for five years. In 2007 the vegetation was cleared and burned again. Information from landholders revealed that the pasture site has been frequently burned over the last 15 years, but the stocking rate has been low with less than one head of cattle per hectare. At the moment of our field campaign, the PA site was covered by grass species Urochloa brizantha and bracken (Pteridium sp.). Our measurements at this site represent the effect of 15 years of intermittent burning and grazing.
The RF site (9°35′48.83″ N, 73°7′2.22″ W) is a secondary old-growth forest that is more than 50-years-old. The RF site was probably burned, similar to the other study sites and there has been selective cutting of trees. Currently, the RF is affected only by the occasional crossing of cattle. We used the RF site as a control area to assess reference values for soil physical and hydraulic attributes.
A graphical summary of the land-use history for the four land-use types described previously is provided in Figure 3.

Vegetation Sampling
In each 79 m 2 plot, we measured the following vegetation attributes: (1) canopy cover; (2) vegetation height; (3) DBH; and (4) total number of trees. These are key ecological indicators, useful to evaluate vegetation structure in tropical forest restoration projects [25,26]. In each plot, we surveyed all living trees with diameter at breast height (DBH) > 5 cm. The vegetation height was estimated with a 3 m measuring stick, and the remaining height of trees taller than this was estimated visually. We furthermore measured the percentage of canopy cover three times with the mobile application CanopyCapture, which is a free, accurate and efficient, canopy cover measurement device [27].

Vegetation Sampling
In each 79 m 2 plot, we measured the following vegetation attributes: (1) canopy cover; (2) vegetation height; (3) DBH; and (4) total number of trees. These are key ecological indicators, useful to evaluate vegetation structure in tropical forest restoration projects [25,26]. In each plot, we surveyed all living trees with diameter at breast height (DBH) > 5 cm. The vegetation height was estimated with a 3 m measuring stick, and the remaining height of trees taller than this was estimated visually. We furthermore measured the percentage of canopy cover three times with the mobile application CanopyCapture, which is a free, accurate and efficient, canopy cover measurement device [27].

Soil Sampling
We selected five sampling points within each study plot for the infiltration measurements. The distance between sampling points was 2 m. For measuring the infiltration rates we used the Beerkan method [28]. In total 40 Beerkan experiments were carried out in February 2021. We used a PVC ring with an inner diameter of 15.5 cm, which was inserted to a depth of about 1 cm into the soil surface. The location of the insertion of the PVC ring was previously prepared by removing the litter and other materials to expose the soil surface. At each measurement point, a known volume of water (150 mL) was repeatedly poured into the ring. Then, the time for each poured volume to infiltrate was logged. This procedure was repeated until the difference in infiltration time between two or three consecutives trials became negligible.
The values of K s (mm h −1 ) were estimated using the Beerkan Estimation of Soil Transfer parameters (BEST) method proposed by Lassabatere et al. [28]. Following Di Prima et al. [29], we used a correction factor to the infiltration model used in BEST to be able to successfully analyze data collected under both hydrophilic (CO, NR, and PA sites) and hydrophobic (RF site) soils. More specifically, it allowed us to account for the rate of attenuation of infiltration rate due to water repellency at the RF site [30].
For three of the five sampling points, we collected undisturbed soil cores (100 cm 3 ) at 0-5 cm depth. For each of these cores, we determined the bulk density, BD (g cm −3 ), and the initial volumetric soil water content θ i (cm 3 cm −3 ). At the end of each infiltration test, a disturbed soil sample was collected to determine the saturated gravimetric water content. The BD was used to calculate the saturated volumetric soil water content, θ s (cm 3 cm −3 ). During the soil sampling we observed soil water repellency phenomenon. We therefore used the water drop penetration time (WDPT) test to quantify the extent of repellency, for the three sampling points at each study plot (six for each land-use type). For each sampling point, the litter and leaf residues were gently removed from the soil surface and ten drops (0.05 mL) of distilled water were placed onto the soil surface using a pipette at a small distance (i.e., a few cm) from each other. We measured the time until complete infiltration of each drop. The recorded time was stopped after 3600 s, although some drops did not infiltrate during this time interval. In such cases, a WDPT value of 3600 s was assigned [31]. Thus, we measured the infiltration of 240 drops, 60 for each site. A representative WDPT value was obtained by averaging the 60 WDPT measurements for each site. The average WDPT values were then grouped in five classes from wettable (WDPT < 5 s) to extremely water-repellent (WDPT > 3600 s) according to the classification proposed by Dekker and Ritsema [32].

Data Analysis
The variables considered in this investigation (K s , BD, θ i , θ s , canopy cover, vegetation height and DBH) were summarized by calculating the arithmetic mean and the associated coefficients of variation (CV) or the standard error for the vegetation variables. For each studied variable, we checked the assumptions of homogeneity of variance based on the Kolmogorov-Smirnov test. One-way analysis of variance (ANOVA) was performed considering the land-use type as an explanatory variable. Then, the Tukey test was applied, which compares the means two by two, the related p-values were computed and compared to the level of significance of 0.05. The K s data were log-transformed prior to analyses, the transformed data were normally distributed according to the Kolmogorov-Smirnov test. Thus, the geometric mean and the associated geometric coefficient of variation were calculated for these variables [33,34]. All analyses were performed in the R software (R Core Team, Vienna, Austria, 2018).

Results
The vegetation characteristics at the CO site differ statistically from NR and RF (Table 1). The canopy cover was significantly different between the CO and NR and RF, highlighting the higher similarity between the NR and RF than CO and NR (Figure 4a). Overall, there was a high similarity between the RF and NR, except for the larger number of trees (32) at NR, than RF and CO (Table 1). Moreover, the average DBH for the RF was smaller than for the younger NR site.

Results
The vegetation characteristics at the CO site differ statistically from NR and RF ( Table  1). The canopy cover was significantly different between the CO and NR and RF, highlighting the higher similarity between the NR and RF than CO and NR (Figure 4a). Overall, there was a high similarity between the RF and NR, except for the larger number of trees (32) at NR, than RF and CO (Table 1). Moreover, the average DBH for the RF was smaller than for the younger NR site.   The results of the canopy cover and some physical soil attributes are represented using box plot graphics (Figure 4). The canopy cover showed significant differences in CO with NR and RF. We highlight the higher similarity between the NR and RF. The K s was significantly lower for the PA (range: 7-20 mm h −1 ) than the other study sites, for which it was similar (range: 71-2127 mm h −1 ) (Figure 4b and Table 2). The BD ranged from 0.72 to 1.42 g cm −3 , there were no significant differences between the plots in terms of BD, although the values were lowest at the RF and the highest at the PA site (Figure 4c). The antecedent moisture content, varied between 0.04 to 0.25 cm 3 cm −3 and was moderately low for PA and significantly different from CO and NR. The θ s varied between 0.15 to 0.63 cm 3 cm −3 with significantly lower average values for NR and RF than the CO and PA sites ( Table 2). Table 2. Minimum (min), maximum (max), mean, and coefficient of variation (CV) of the saturated soil hydraulic conductivity, K s (mm h −1 ), dry bulk density, BD (g cm −3 ), initial volumetric soil water content, θ i (cm 3 cm −3 ), and saturated volumetric soil water content, θ s (cm 3 cm −3 ). CO, shade-grown coffee; NR, Natural regeneration; PA, pasture; RF, Reference forest. Note. For each soil indicator different superscript letters denote statistically significant differences between land-use types (p < 0.05).
The soil water repellency varied among the different land-use types as well and ranged. The results of the WDPT show the existence of different degrees of soil water repellency ranging from wettable to extremely water repellent. The RF was the only landuse with severe and extreme water repellency, 22% and 10% of all the drops, respectively. Strong water repellence was detected at RF (60%) and NR (3%), and the slightly water repellent-class was observed in all land-use types. The wettable class was dominant at PA (98%) and CO (87%) sites ( Figure 5). with NR and RF. We highlight the higher similarity between the NR and RF. The Ks was significantly lower for the PA (range: 7-20 mm h −1 ) than the other study sites, for which it was similar (range: 71-2127 mm h −1 ) (Figure 4b and Table 2). The BD ranged from 0.72 to 1.42 g cm −3 , there were no significant differences between the plots in terms of BD, although the values were lowest at the RF and the highest at the PA site (Figure 4c). The antecedent moisture content, varied between 0.04 to 0.25 cm 3 cm -3 and was moderately low for PA and significantly different from CO and NR. The θs varied between 0.15 to 0.63 cm 3 cm -3 with significantly lower average values for NR and RF than the CO and PA sites ( Table 2). Note. For each soil indicator different superscript letters denote statistically significant differences between land-use types (p < 0.05).
The soil water repellency varied among the different land-use types as well and ranged. The results of the WDPT show the existence of different degrees of soil water repellency ranging from wettable to extremely water repellent. The RF was the only landuse with severe and extreme water repellency, 22% and 10% of all the drops, respectively. Strong water repellence was detected at RF (60%) and NR (3%), and the slightly water repellent-class was observed in all land-use types. The wettable class was dominant at PA (98%) and CO (87%) sites ( Figure 5).

Discussion
Forest restoration and agroforestry in tropical environments of South America are receiving attention for their ecological, economic and social implications [8,35]. Our investigation, focused on the study of the soil hydraulic properties in Sub-Andean Forest. It revealed that the K s at the CO plantation was statistically similar to that of NR and RF. The lower K s at NR did not allow us to fully accept our study hypothesis, namely that averages of the K s values would vary as follows: RF > NR > CO > PA. Instead, the order is RF = NR = CO > PA, despite the significantly lower canopy cover at CO site. The high K s at CO plantation can be possibly due to litter thickness [36] and extensive root systems created by the woody coffee plants [37]. Coffee agroecosystems in the topsoil have a greater abundance of medium to fine roots, while the rhizosphere of secondary old-growth forests is dominated by the thick roots of mature trees and greater abundances of medium roots [36,38]. Another important factor that helps to explain the high K s at CO plantation are the management activities. For example, at our CO site the coffee bushes are pruned regularly, leaving the residues on the topsoil, in according to other studies this practise can favor the hydrophysical soil attributes [37,39]. Moreover, our CO plantation has a low-intensity use, which is characterized by the absence of tillage and soil compaction, the few soil disturbances at CO site are associated mainly with human trampling. In contrast when the intensity of the land-use increases soil infiltration can be negatively affected. Cannavo et al. [40] on Costa Rican Andisols found that the soil infiltration was lower at coffee monoculture with high-intensity use than for shade-grown coffee. These findings emphasize the need to experimentally check the real soil conditions that arise, over time, in agroforestry environments.
The canopy cover at CO site was significantly lower when compared with NR and RF. The low canopy cover in CO can be associated to the few tree individuals (i.e., 5 trees) and bananas (Musa spp.). The shade tree species in the agroecosystems can play an important role for soil water infiltration [41,42]. In particular, leguminous trees (Erythrina spp. and Inga spp.) that are commonly used in coffee plantations can increase soil water infiltration [37,40]. We suggest that more research is necessary to understand tree-species effects on hydrological and environmental services in the coffee agroecosystems.
The shade trees in coffee agroecosystems can provide important benefits such as: (i) decrease the loss of topsoil by reducing erosion [4,37], (ii) limit water loss during the dry season [43,44], (iii) provide nitrogen via litter decomposition [36], (iv) regulate pests and diseases [45], (v) increase carbon sequestration [46], (vi) provide a better quality of coffee beans, and (vii) favor biotic diversity [36,47]. On the other hand, a high canopy cover in agroforestry coffee plantations can affect the yield, representing financial loss for coffee farmers [40]. Thus, the pros and cons of shade-grown coffee in our study area need to be evaluated in future studies, collecting more wide-ranging empirical data.
As expected, the mean K s value was lower in PA (12 mm h −1 ) than in forests stands (244-469 mm h −1 ). This result can be related to the highest BD found in the pasture plots. Also, our result agrees with several other studies [48][49][50]. In tropical pastures, some of the factors that affect K s negatively are the soil compaction, the intensity of the cattle-grazing, the lower soil faunal activity and organic matter [10,48,51].
Assuming that the K s of degraded land was similar to the measurements at the PA site prior to the regrowth of the forest, our results show that abandonment and forest regrowth on degraded land can improve the K s , and that after 15 years the K s became similar to that of the reference forest. The increase of K s during the natural regeneration process has been reported by several author before [10,26,[51][52][53]. Ziegler et al. [49] suggest that a period of 25 years could be enough to recover K s after 2-4 years of swidden cultivation in Vietnam. Leite et al. [52] mention that after grazing, the complete recover of the K s will probably take more than 35 years in the Caatinga ecosystem. Zwartendijk et al. [53] point out that the full hydrological recovery of tropical degraded sites with a natural regeneration process, may take several decades. However, K s does not always increase during natural regeneration, in some cases the evolution of K s and hydrological functions can be a relatively slow process [15,50,54,55]. Pereira et al. [15] showed for the Atlantic forest biome that after 42 years the natural regenerated forest still had similar K s values as the sugarcane site [15]. They attributed the results to the legacy effect of intense land degradation before forest regrowth. Similarly, Ghimire et al. [55] reported for a 25-year-old forest in Nepal similar K s values as for pasture sites. They explained that removal of litter, grazing and fuelwood harvesting affected K s recovery. In summary, the soil recovery process after natural regeneration or tree planting depends on many factors and is still not fully understood. The multi-decadal effects of vegetation regrowth has not been captured in most of the studies [56,57].
At our sites, the highest K s values were found at the RF, which was also the site with the lowest BD and greater values for the vegetation attributes (i.e., canopy cover, vegetation height of tree, DBH and total number of trees). This result can be explained by the more conserved soil conditions and the longer recovery process of the vegetation (i.e., more than 50 years). Many other studies have observed a high K s for little disturbed tropical forests [53,55,56,58]. Our result thus supports the importance of preserving the forest cover in transformed landscapes in order to maintain infiltration process, soil hydraulic efficiency and, therefore, the regulation of the hydrological cycle. In addition, the K s values obtained in the RF plots reflect the high spatial variability of the infiltration. Although we did not directly quantify variables such as soil macroporosity, soil biological activity, root biomass, plant diversity, organic matter and topographic variations are generally considered the main drivers of the spatial variability of K s in forest soils [26,52,53,55,56]. Given the importance of the aforementioned factors, they can be the subject of specific investigations in the future.
Soil water repellency is a transient soil property poorly understood for tropical soils [59,60]. The WDPT allowed to detect a stronger water repellency at the NR and RF sites. Our findings are in agreement with previous studies in tropical ecosystems [61][62][63], which have reported hydrophobicity under forest cover. It is important to highlight that soil water repellency can negatively affect the soil infiltration process and increase soil erosion [48]. Soil water repellency can be attributed to many factors such as amphiphilic molecules produced by plants and organisms, soil conditions after forest fires or dry periods. Soil texture, soil temperature, pH, soil water content and soil organic carbon can affect soil water repellency as well [59,60]. Future research needs to improve the understanding of factors promoting the soil water repellency in tropical forest soils. In Colombia, the soil water repellency has been reported for Andisols with Pinus patula plantation [64], our work reports the first measurements for a Sub-Andean Forest remnant.

Conclusions
From this study it can be concluded that agroforestry systems such as the shade-grown coffee in Sub-Andean Forests benefit the hydrophysical soil attributes. Our study also concludes that abandonment and forest regrowth on degraded land can improve the K s , and that after 15 years the K s became similar to that of the reference forest. The comparative analysis of K s in agroforestry use (CO), secondary forests (NR and RF) and pasture (PA) show the similarity in K s between the land-uses with trees (CO, NR, and RF). Furthermore, the PA site with a higher intense land-use history showed the lowest K s. Funding: This research was funded by the Colombian Ministry of Science (grant no. 848) and "Vicerrectoría de Investigaciones, Extensión y Proyección Social" of Universidad del Atlántico.