Hydrological Effects of Spatial Harvest Patterns in a Small Catchment Covered by Fast-Growing Plantations in the Neotropics

Lara Gabrielle Garcia ^1,*, Luiz Felippe Salemi ², Rafael Mingoti ³, Carla Cristina Cassiano ¹, Aline Aparecida Fransozi ¹, Vinicius Guidotti de Faria ⁴, Carlos Alberto Vettorazzi ⁵, Walter Paula Lima ¹, Silvio Frosini de Barros Ferraz ¹

1. Forest Science Department, University of São Paulo/ESALQ, Av. Padua Dias, 11, 13400-900, Piracicaba, São Paulo, Brazil

2. Planaltina Faculty, University of Brasília (UnB), Vila Nossa Senhora de Fátima, Planaltina, Brasília-DF, 73345-010, Brazil

3. Territorial Management, Brazilian Agricultural Research Corporation (Embrapa), Av. Soldado Passarinho, 303, 13070-115, Campinas, SP, Brazil

4. Institute of Forest and Agricultural Management and Certification (IMAFLORA), Estrada Chico Mendes, 185, 13426-420, Piracicaba, São Paulo, Brazil

5. Biosystem Engineering Department, University of São Paulo/ESALQ, Av. Padua Dias, 11, 13400-900, Piracicaba, São Paulo, Brazil

* Correspondence: Lara Gabrielle Garcia

Academic Editor: Wen-Cheng Liu

Special Issue: Advances in Hydrology, Water Quality and Sediment Simulation Modelling

Received: October 07, 2022 | Accepted: February 16, 2023 | Published: February 21, 2023

Adv Environ Eng Res 2023, Volume 4, Issue 1, doi:10.21926/aeer.2301022

Recommended citation: Garcia LG, Salemi LF, Mingoti R, Cassiano CC, Fransozi AA, de Faria VG, Vettorazzi CA, Lima WP, de Barros Ferraz SF. Hydrological Effects of Spatial Harvest Patterns in a Small Catchment Covered by Fast-Growing Plantations in the Neotropics. Adv Environ Eng Res 2023; 4(1): 022; doi:10.21926/aeer.2301022.

© 2023 by the authors. This is an open access article distributed under the conditions of the Creative Commons by Attribution License, which permits unrestricted use, distribution, and reproduction in any medium or format, provided the original work is correctly cited.

1. Introduction

Eucalyptus forests are associated with potential impacts on water resources due to their high water use [1,2,3,4,5]. For this reason, countries like Brazil with approximately 7.5 million hectares [6] need to understand the effect of such forest plantations.

Existing research recognizes the critical role forest management practices play in tackling their impact on water resources. For example, using forest stands of different ages [7,8] and increases in rotation periods [9] are among the management action taken to reduce the hydrological effect of Eucalyptus plantations [10]. Another way to cope with this problem is to manage two things: (i) the amount of harvest area and (ii) the hillslope position where the forest is being harvested. Regarding the former, various papers have already shown the positive relationship between the amount of forest harvested area and streamflow (e.g. [1,2]). As for the latter, due to the differences in hydrological processes along topographic gradients [11,12,13,14] and the main impact of forest harvest in hydrological processes [2,15,16,17], forest harvest position could be used as a management strategy for mitigating forest plantation negative effects on streamflow. For example, owing to differences in evapotranspiration along the hillslope [18], previous simulations have shown that forest harvest carried out near stream areas could generate higher streamflow than the same amount of harvest area carried at the catchment divide [11]. However, previous studies focused on natural Eucalyptus forests [18] and temperate forests [11]. Thus, important information regarding Eucalyptus plantations which, due to their high growth rates, could have more pronounced effects than other forest types, is still lacking. Thus, understanding such effects within Eucalyptus plantations could add more knowledge on the hydrological effects of forest management and, consequently, provide more insights to formulate management actions that could minimize their hydrological effect.

In this context, the present paper had a two-fold objective: (i) assess the effect of spatial harvest position (downslope × upslope position) on streamflow and (ii) understand the combined effect of rainfall amounts and forest harvest on peak and volume flows.

2. Materials and Methods

2.1 Hydrological Data

The study area is an experimental catchment located in Itatinga, Sao Paulo State (23°02’01’’ S; 48°37’30’’ W, altitude of 850 m) (Figure 1). The region is representative of the climatic conditions and management practices in fast-growing plantation forests in Brazil. The catchment is 85.8 ha and is covered by Eucalyptus and Pinus forests (85%), riparian native vegetation around the stream channel (10%), and roads (5%).

Figure 1 Location of study area: São Paulo State in Brazil and Tinga catchment land use (23°02’01’’ S; 48°37’30’’ W, altitude of 850 m).

The mean annual temperature is 19.9°C, with the mean annual precipitation ranging from 850 mm to 1600 mm, which is most concentrated in the summer (October to March), with an average water surplus of 389.8 mm and a deficit of 2.4 mm [19]. The soil types present in the study area are Typic Hapludox and Rhodic Hapludox [20], and the soil texture varies from sandy clay to sandy clay loam [20].

Streamflow data were collected from a gauge using an H type flume equipped with an automatic stage recorder (Thalimedes Shaft Encoder sensor, OTT), with a 15-minute resolution, coupled with a data logger. Precipitation data were measured with an automatic rain gauge (TR-525I, Texas Electronics, tipping bucket, 0.5 mm resolution) located 1 km from the stream gauge. This gauge is coupled with a datalogger and records rainfall amounts at 30-min intervals. Other meteorological data used include net radiation (MJ m^-2 d^-1), air temperature (°C), air humidity (%) and wind speed at 2 m (m s^-1). These data were obtained from an eddy flux tower located 13 km from the catchment at 30-min recording intervals.

2.2 Hydrological Modeling

The methodological approach used in this study simulates forest harvest scenarios using a physical-based distributed hydrological model (Gridded Surface Subsurface Hydrologic Analysis - GSSHA) [21] to understand the effects of (1) forest harvest amount (70% and 30% of the total catchment area) on water yield and flow regime; (2) harvest spatial position (downslope × upslope position) on water yield and flow regime; and (3) combined effects of rainfall amounts and forest harvest on peak and volume flows. This study collected data from an experimental catchment in Sao Paulo state, Brazil.

The GSSHA code was used to carry out the hydrological modeling in the catchment to simulate the hydrological process following land use/land cover changes. Both forest plantation effects and forest cover changes have been modeled with several hydrological models [22,23,24], including distributed models under tropical conditions [25]. However, to our knowledge, this is the first use of the GSSHA in tropical conditions and with a fast-growing Eucalyptus forest.

The GSSHA code is a spatially distributed, process-based model and has been developed to simulate processes in hydrological cycles at a catchment scale in continuous simulations or rainfall events [21,26]. The model has been used to simulate runoff production mechanisms [21,27], land-use pattern changes [28], impacts of land-use spatial distribution changes [29], and water flows on forest catchments [30].

2.3 Model Parameterization

The land-use storage capacity was derived from the relation between the leaf area index (LAI m²·m^-2) and the interception index (α), assumed to be 0.2 mm [31]. LAI values for Eucalyptus and Pinus plantations were obtained from Almeida et al. [32] and Cabral et al. [33]. The channel cross-section was defined assuming a trapezoidal area with 2:1 (H:V) side slopes, 0.5 m bottom widths, 0.5 m depth and a roughness coefficient of 0.035, and a 1D flow model was used for channel routing. Parameter values are presented in Table 1, and model calculations were performed at a time step of 10 seconds.

Table 1 Parameter values used in GSSHA model for two soil type and five land uses in Tinga catchment.

The scenarios were run using a uniform grid (spatial discretization: 5 × 5 m resolution) created from a catchment digital elevation model (5 m resolution). The elevation model was used to compute the flow direction, generate the channel network and delineate catchment boundaries.

Since long-term simulations were used, the infiltration process was computed using the Green and Ampt with redistribution (GAR) method (1D) [44]. This method is particularly useful in studying fine soil textures such as those in this study catchment [20] since it is a method that is in good agreement with the Richard equation and allows the simulation of saturation excess overland flow [44].

Soil hydraulic properties were obtained by laboratory methods, using soil samples collected in the field. To determine the hydraulic conductivity and soil water retention curves, a constant head permeameter, and pressure plate extractor were used, respectively (Table 1).

Overland flow routing was computed using the Alternating Direction Explicit (2D overland flow) algorithm that requires Manning’s roughness coefficient (n) as an input parameter for the four land uses. We derived values of n from the Curve number coefficient (CN) through Eq. (1) [45] and Eq. (2):

\[ T_c=\frac{l^{0.8}\left[\left(\frac{1000}{CN}-10\right)+1\right]^{0.7}}{1140Y^{0.5}} \tag{1} \]

\[ n=\frac{\left(\frac{T_cP_2^{0.5}S^{0.4}}{0.007}\right)^{\frac{1}{8}}}{l} \tag{2} \]

where

T_c = time of concentration, h

l = flow length, ft

CN = curve number factor

Y = average watershed land slope, %

n = Manning’s roughness coefficient

P₂ = 2-year, 24-hour rainfall, in

S = slope of land surface, ft/ft

For this method, we must take method limitations into account [45] and consider flow length (l), equal sheet flow length (l_S), average watershed land slope (Y), and equal land surface slope (S) multiplied by 100 (S*100); P₂ equals 1.651 inches. CN for Eucalyptus, Pinus, riparian native vegetation and roads are, respectively, 25, 45, 20 and 72 [42,46]. The n values used in the model for Eucalyptus, Pinus, riparian native vegetation and roads are, respectively, 0.132, 0.149, 0.346 and 0.046 m s^-1/3.

The potential evapotranspiration was calculated using the Penman-Monteith method. Parameters like land use, albedo, vegetation transmission coefficient, canopy stomatal resistance, and vegetation height were obtained from previous studies (see Table 1). Eucalyptus height was derived from an allometric equation in Christina et al. [47]. Actual evapotranspiration is calculated as a function of the vegetation cover and its root system and canopy resistance characteristics [21].

Lateral groundwater flow was simulated using a 2D groundwater flow model (Trescott-Larson method) with groundwater and stream channel following Darcy’s law [21].

2.4 Model Calibration and Validation

The simulations used warm-up conditions (from January 2010 to December 2010). The model was calibrated for a continuous period from January 2011 to April 2012, when the catchment was entirely covered by forest (native and planted). The model was manually calibrated against the daily observed flow for a continuous simulation period.

The model was calibrated and validated to maximize the modified Kling-Gupta efficiency index (KGE’) and its decomposition [48,49], the Nash-Sutcliffe efficiency (NSE) [50] and the efficiency index of Bravais-Pearson (R²) [51]. The modified Kling-Gupta efficiency index (KGE’) was decomposed into three components (correlation coefficient r, variation coefficient ratio γ_KGE’ and relative bias β_KGE’) [49]. A value of 1 for all coefficients is considered the best agreement between the observed and predicted data.

We used resampling techniques as the K-fold cross-validation method to cross-validate the model [52]. The flow period data were split into three periods where two of them were used for calibration, and the remaining period was used for validation (split 1: Jan/2012 to Mar/2012; split 2: Apr/2012 and Jan/20122; split 3: Feb/2011 to Apr/2011). For an additional check, we used flow duration curves (FDCs) to assess the differences in flow duration curves and total flow between the predicted and observed periods as a supplementary metric for model calibration [49,52].

2.5 Harvesting Simulation

First, forest harvest amounts were simulated using scenarios of 30% and 70% of the harvest area (hereafter H30 and H70, respectively) of Eucalyptus and Pinus plantations. Therefore, native vegetation present in the catchment was not included. The scenarios used to represent typical forest management procedures in Brazil [7].

Second, to clarify the effect of the harvest position, each scenario (H30 and H70) was combined with spatial patterns (downslope and upslope). Those scenarios were referred to as H30_DOWN, H30_UP, H70_DOWN, and H70_UP (Figure 2). Land use after harvesting was considered the “harvest area” and corresponded to bare soil.

Figure 2 Harvest scenarios used to simulate amount of harvested areas and spatial patterns: (a) reference scenario with 100% forest; (b) 30% harvest at a downslope position (H30DOWN); (c) 30% harvest at an upslope position (H30UP); (d) 70% harvest at a downslope position (H70DOWN) and (e) 70% harvest at an upslope position (H70UP).

To analyze forest harvest amounts and spatial patterns, three flow indices, namely, Q70 (low flow), Q50 (median flow), and Q10 (high flow) and the annual water yield obtained from FDC, were considered.

Finally, the effect of the spatial pattern on peak flow was tested under two rainfall intensities. The first intensity, R5, was defined as the 5^th percentile of the highest recorded rainfall during the predicted period (P = 50.5 mm during 5 hours with a maximum of 30 minutes of intensity at 59.4 mm h^-1). In the same way, the other rainfall, namely, R95, is the 95^th percentile of the lowest rainfall recorded in the area during the predicted period (P = 8.1 mm during 3.5 hours with a maximum of 30 minutes of intensity at 21.8 mm h^-1).

Data analysis was performed by comparing the results of the harvest scenarios with a reference scenario of 100% forest cover (i.e., 0% harvest). Hydrological simulations were carried out for one year (72 rain events). The first three months were considered a warm-up period and were not included in data analyses. The parameters modified to these simulations were: Land Surface Albedo (0.25); Vegetation Height (0.0 m); Vegetation transmission coefficient (1.0); Canopy Stomatal Resistance (545.2 s m^-1); Storage capacity (0.0 mm); Interception coefficient (0.0 mm h^-1); Retention depth (1.0 mm); CN (72).

3. Results

3.1 Calibration and Validation

The hydraulic conductivity of the saturated zone appeared to be the most sensitive parameter in the model [21]; therefore, it was used in the calibration processes. After calibration, the values of this parameter were 0.29 and 1.2 cm h-1 for the sandy clay and sandy clay loam soils, respectively.

The calibration result is given in Table 2. The coefficients index NSE and KGE’ show an agreement of 0.6 for predicted and observed daily data, which was satisfactory goodness of fit [53]. Perhaps, the values of NSE at approximately 0.6 indicate that, for some periods, the model could not predict the observed values with the good agreement [52]. The bias coefficients indicate a 5% overestimation of the predicted mean flow data about the observed flow data. The predicted flow had a linear correlation of 80% with the observed flow.

Table 2 Results for the calibration and validation using the K-fold cross-validation method. The numbers 1, 2 and 3 refer to the split data thirds (split 1: Jan/2012 to Mar/2012; split 2: Apr/2012 and Jan/20122; split 3: Feb/2011 to Apr/2011).

Flow duration curves showed a better agreement for the high and mean flows than for the low flows (Figure 3). However, total flows for the entire period (141 mm and 149 mm for observed and predicted values, respectively) had lower than 5% between predicted and observed values.

Figure 3 Predicted and observed flow duration curves for daily flow data.

3.2 Effect of Spatial Harvest Patterns on Water Yield and Flow Regime

Harvesting spatial patterns could influence flow regimes since forest harvesting on downslope positions (in both H30_DOWN and H70_DOWN scenarios) presented greater FDC indices and water yield than those at upslope positions (Table 3).

Table 3 Summary of the effects of the harvest scenarios at downslope (H30_DOWN and H70_DOWN) and upslope (H30_UP and H70_UP) positions.

Regarding effects on water yield, the H30_DOWN scenario resulted in a 12% increase relative to the same harvest area on the upslope (H30_UP). The same trend was found for H70 scenarios, that is, the effect of the H70_DOWN on water yield was 11% greater compared to H70_UP.

Flow regime indices showed that the low flow index increased more than 17% for the H70_DOWN, whereas the H70_UP yielded only 4%, and the same trend was found for medium and high flows. Regarding spatial harvest patterns, the H30_DOWN scenario promoted a greater increase in FDC indices (increases of 14%, 12% and 6% for low, medium and high flow, respectively), when compared to the H70_UP scenario (increases of 4%, 4% and 3% for low, medium and high flow, respectively) (Figure 4). The H30_UP scenario FDC indices did not differ from the reference scenario (Table 3).

Figure 4 Flow duration curve indices from both H70 and H30 scenarios. Reported effects are relative to the reference scenario (100% forest cover).

3.3 Combined Rainfall and Forest Harvest Effect on Peak Flow

Peak and volume flows were affected for all harvest amount scenarios since harvest increases peak and volume flows at R5 and R95 events. However, the results under R5 presented greater absolute increases in peak flows compared to R95 for both the H30 and H70 (Figure 5 and Figure 6).

Figure 5 Peak flow for the two rainfall amounts (R5 and R95) in the 30% harvest scenarios down (H30_DOWN) and upslope (H30_UP) and in the 70% harvest scenarios down (H70_DOWN) and upslope (H70_UP). H0: Reference scenario (0% harvest).

Figure 6 Volume flow for the two rainfall amounts (R5 and R95) in the 30% harvest scenarios down (H30_DOWN) and upslope (H30_UP) and in the 70% harvest scenarios down (H70_DOWN) and upslope (H70_UP). H0: Reference scenario (0% harvest).

For the H30 under R5 and R95, the increase in peak flows is greater for the H30_DOWN than for the H30_UP, but that was not the case for the H70, where peak flow under the R5 showed lower values for downslope compared to the upslope scenario (Figure 5), resulting in a small difference of 3% (or 1.1 L s^-1).

Under R95, the H70_DOWN and H70_UP had greater effects on increasing peak flows (21% and 13% relative to the reference scenario, respectively) compared to the H30 in the same positions (18% and 1% relative to the reference scenario, respectively) (Figure 5 and Figure 6). Nevertheless, for R5, our results showed that the H30_DOWN generated a greater increase in peak flow relative to the H30_UP (Figure 5 and Figure 6). In addition, for the H70, small differences were detected between the upslope and downslope harvests.

4. Discussion

4.1 Forest Harvest Spatial Patterns: Effect on Water Yield and Flow Regime

We showed that spatial harvest patterns could influence water yield, with a greater increase in water yield after a forest harvest at downslope compared to upslope positions. These results align with similar studies elsewhere [11,54]. For example, Abdelnour et al. [11] found that a 20% harvest of conifers near a stream (downslope position) resulted in a water yield increase of 73% relative to a harvest on the upslope position. According to the authors, these results may be because subsurface flow generated from the upslope position, as opposed to the downslope position, has a longer flow path to reach the stream channel. This explanation seems to be related to the increased time that water is available for trees as they move downslope. Thus, one of the explanations for the increase in water yield could be the changes in the other hydrological processes (interception, transpiration and overland) [2,55].

Previous studies have shown that vegetation at downslope may have higher water use than vegetation at upslope positions [55,56]. These differences can be explained, in part, by the proximity to the water table in the downslope position, which could support the hypothesis that evapotranspiration next to streams is greater [57]. Upslope Eucalyptus natural forests had 40% lower evapotranspiration rates than the middle and downslope [18]. However, the study by Mitchell et al. [18], may be somewhat limited by variations in soil, vegetation type (natural Eucalyptus forests), and density across the hillslope; thus, the effect of slope positions were obscured by these other factors. In this study, we modeled the same soil and vegetation across the hillslope. Therefore, the effect of the hillslope position can be considered dissociated from other influencing factors. Yet, the results are in line to some degree.

As for flow regime, our results indicate that effects on FDC indices (low, medium and high flows) follow the same trend observed for water yields. Forest harvests and their spatial locations can also modify flow regimes. These findings suggest that spatial patterns of harvests could influence flow distribution; therefore, forest harvest should be considered in decisions to improve water availability. Forest plantations will generate lower water yield reductions when located in areas far from the stream [56] due to lower evapotranspiration rates [18].

4.2 Combined Rainfall and Forest Harvest Effects on Peak and Volume Flow

Some authors report that the influence of forest cover on streamflow generation depends on rainfall event characteristics, which is relevant for high rainfall intensities and amounts [54,58,59].

Our results showed that all harvest scenarios (different amounts and spatial positions) increased peak and volume flow for both rainfall amounts (R5 and R95 events). As expected, low rainfall amounts (e.g., R95) had lower increases than high rainfall amounts (R5). These results are by previous studies [54,60,61].

Under low rainfall amounts, a 70% harvest results in greater peak and volume flow than 30% harvest scenarios. Although the harvest amount effects on peak and volume flows were more noticeable, harvests occurring in the downslope positions still showed a larger effect than the same harvest amount scenario in the upslope counterparts. Similar results have been found elsewhere (see Abdelnour et al., [11]).

For high rainfall events (R5), 70% of forest harvest scenarios did not differ in peak flows between the upslope and downslope (an increase of 19% about the reference scenario for both spatial positions). The likely explanation for the such finding is that, under extreme events, harvest amount has little effect on peak flows [61,62]. For the 30% forest harvest, the scenarios resulted in an increase of 18% and 6% for the downslope and upslope, respectively. The reason for such a result might be related to the longer path water has to take to reach the stream when the harvest is carried out in upslope areas. Upslope water moving downward may be subject to evapotranspiration to a greater degree compared to downslope water.

Differences in peak flow for harvest amount scenarios (30% vs 70%) could be observed when harvest occurs at upslope (6% and 18% increase in volume flows relative to the reference scenario, respectively). Once again, these results are likely related to harvest-suppressed transpiration rates that increase the water stored in the soil [61] readily generating quicker flow. As more water is available to flow at 70% harvest, more water runs downward generating greater volume flows.

4.2.1 Management Implications

Based on our results, we suggest two strategies to minimize the effect of Eucalyptus fast-growing forest plantations on water resources:

1. (i) Where water availability is lower (due to a drought event or a normal lower water availability region), harvest amount should always be greater compared to more humid regions, generating more water available for downstream users.
2. (ii) Regarding harvest position, forest harvest in downslope areas is also indicated in dry regions since more water will be available to streamflow.

5. Conclusions

Harvest spatial patterns had different influences on flow metrics where harvest at downslope had greater increases in water yield compared to the same harvest at upslope positions.

In terms of sustainability in Eucalyptus plantations, the incorporation of forest harvest amounts and spatial positions into management strategies could be a useful tool to mitigate hydrological effects and reconcile biomass production and water resources conservation. These results directly affect management strategies, particularly in water-limited regions, and provide insight into how to mitigate Eucalyptus plantations’ negative effects on water resources.

Acknowledgments

This work was supported by FAPESP grants n°2012/07428-3, 2016/02877-5. We gratefully acknowledge the team of Itatinga Forest Science Experimental Station and the team of Forest Hydrology Lab ESALQ/USP.

Author Contributions

Lara Gabrielle Garcia – conceived and designed the analysis, collected the data, performed the analysis, wrote the paper. Luiz Felippe Salemi - conceived and designed the analysis, wrote the paper. Rafael Mingoti - conceived and designed the analysis. Carla Cristina Cassiano - conceived and designed the analysis, collected the data. Aline Aparecida Fransozi - conceived and designed the analysis, performed the analysis. Vinicius Guidotti de Faria - conceived and designed the analysis, performed the analysis. Carlos Alberto Vettorazzi - conceived and designed the analysis. Walter Paula Lima - conceived and designed the analysis, collected the data. Silvio Frosini de Barros Ferraz - conceived and designed the analysis, collected the data, funding.

Competing Interests

The authors declare that they have no conflict of interest.