Introduction 

The relationship of vegetation cover and groundwater resources has drawn considerable scientific attention over the last decades. Many studies have shown that deforestation by logging or of natural origin (forest fire, wind damage) increased the average runoff from the affected area (e.g., [4]) and afforestation decreased runoff. Similar results were detected in recharge rates and groundwater depth, predominantly in dry regions. In Australia changes of vegetation from woodland to grassland or crops resulted in increases in recharge rates of one to two orders of magnitude ([22]). The changes were largely due to the altered interception loss and the increased water extraction from the root zone. However, results of paired watershed research agreed only on the direction of the changes, not on their magnitude ([2]).

Comparative water balance studies of forest and low vegetation covers have generally shown higher water use of forest cover ([18], [27], [33]). Nachabe et al. ([27]) analyzed the groundwater consumption in a shallow water table environment and estimated the annual ET for a forested area (1320 mm) and for a pasture (700 mm) using detailed soil moisture and water table monitoring. At the same time, a few studies have found negligible differences in evapotranspiration of different vegetation covers (e.g., [32]). The latter study found that, although there were seasonal differences, on an annual basis the drainage below broadleaved woodland did not significantly differ from that below a pasture.

Due to climate change, air temperature is expected to rise significantly during this century ([14]). As a direct consequence of warmer temperatures, the hydrological cycle will undergo changes with accompanying alteration in the rates of precipitation and evaporation. In Hungary, summer temperatures is expected to increase by up to four degrees by the end of this century, while precipitation is likely to increase in winter and decrease in summer ([10]).

The impact of climate change on groundwater resources was reviewed lately by Green et al. ([11]). Although the uncertainty of predictions for change in groundwater recharge rates and discharge is large (e.g., future climate scenarios and groundwater extraction), numerous studies found that groundwater resources appear to be threatened by future climate change due to increased natural and human water demand.

In light of the water balance uncertainties and the increasing pressures on groundwater resources due to future climate change, a comparative water balance study of an oak forest and fallow vegetation plots was initiated in a drought-threatened lowland environment in north-eastern Hungary. Water balance components were estimated by the Hydrus 1-D numerical model ([34]), calibrated on measured soil moisture and groundwater levels.

Our research questions were the following: (1) what is the magnitude of evapotranspiration components in the two different land-uses; and (2) how has groundwater consumption evolved during dry and wet growing seasons?

  Materials and methods 

Study area

The study plots are situated in the north-east part of Hungary at latitude 47° 58’ N and longitude 21° 42’ E (Fig. 1), built up from sandy river deposits of the early Pleistocene ([3]). The climate is continental; the mean annual precipitation (years 1951-2000) is approximately 520 mm and mean monthly temperature (1951-2000) ranges between -2.4 °C in January and 20.5 °C in July.

The plots were located on the discharge area of the local phreatic groundwater flow system with shallow groundwater levels and groundwater supplies from the adjacent areas. Both study plots were selected so that the elevation of the plots above the nearby ditch was almost the same. The generally similar site conditions of the field plots made it possible to compare water balance components and groundwater consumption. The approximately 300 m distance from the nearby ditch minimized any effect that floods had on water table levels. Surface runoff was not observed since the slope of the plots was less than 2‰ (Tab. 1). For detailed description of the oak plot, see Móricz ([26]).

The naturally regenerated oak forest had 60% pedunculate oak (Quercus robur) which are 20-25 m high and the density is 270 trees ha^-1. The fallow plot, situated about 3 km from the oak plot, is part of an agricultural field that had been plowed in the past, but it has not been seeded and cultivated now for several years and is currently under natural succession. The vegetation consists of furrow-weed with species like Solidago gigantea, Artemisia vulgaris, Amaranthus retroflexus and Ambrosia artemisiifolia.

The vertical distribution of the root system was surveyed in situ by taking three replicate volumetric soil core samples at six depths (0-0.2 m, 0.2-0.4 m, 0.4-0.6 m, 0.6-0.8 m, 0.8-1.0 m and 1.0-1.2 m) for both study plots. After separating the fine roots (diameter < 2 mm) by sieving the soil samples, they were scanned and the total root extent of each sample was determined by grid-counting. Finally, the fine root fraction of each layer was related to the total root extent of the profile. Fine roots reduced approximately linearly with depth at the oak plot (estimated root depth: 1.5 m) and logarithmically at the fallow plot (estimated root depth: 0.8 m).

Soil analyzes included the sieving and hydrometer analyzes of particle size distribution of the soil samples, taken at 0.2 m intervals down to a depth of three meters. At both plots the soil texture was compacted fine sand (0.02-0.2 mm) close to the surface varying between 80-99%. The clay and silt fractions were high below one meter depth at both plots and reach 30-35% and 20-25%, respectively.

Three repetition of undisturbed soil samples for water retention were analyzed at depths of 0.1, 0.3, 0.5, 0.7, 0.9 and 1.2 m by cylinders of 100 cm³ (Tab. 2).

The maximum of the Leaf Area Index (LAI - m² m^-2) was estimated by collecting leaf litter on the ground. At the oak plot, newly defoliated leaves were collected carefully from five locations (1 x 1 m) during late autumn of 2007 and dried in an oven (105 °C for 24 hours). Determination of LAI included the calculation of the ratio of weight to leaf area for a subset of leaves and then for the whole sample. The average LAI of the samples was 3.9 m² m^-2. The 16-day Enhanced Vegetation Index (EVI) product of MODIS (Moderate Resolution Imaging Spectroradiometer) was used to describe the seasonal change of LAI. The 250 m resolution EVI was converted to LAI using the relationship proposed by Wang et al. ([42]).

At the fallow plot, the maximum leaf area index was estimated by leaf collection from three locations (0.5 x 0.5 m). All the collected leaves were scanned and the leaf area was determined by grid counting. The mean LAI of the three samples was 1.1 m² m^-2. The same leaf area index was used throughout the vegetation period. During the dormant season, we assumed the LAI 0.5 m² m^-2, based on the biophysical parameter table of [35] - Fig. 2).

The albedo was derived from the 16-day estimates of the MODIS images. The missing values of albedo in winter were assumed as snow cover and were replaced by an albedo of 45% at the oak and 75% at the fallow plot ([17]).

Monitoring at the plots

Meteorological variables, soil water content and groundwater level were monitored by automatic equipment.

A weather station (iMETOS, Pessl Instruments, Austria) at the fallow plot and an automatic rain-gauge (Rainlog Data Logger, Rainwise, USA) 500 m apart from the oak plot were employed to monitor weather conditions (Fig. 1). The time resolution of rain measurements was 1 minute, while all other variables were recorded hourly.

The volumetric water content was monitored with FDR (Frequency Domain Reflectometers) using Decagon EC-5 probes (Decagon Devices, Pullman, USA) with a time interval of 15 min. The probes were installed at depths of 0.1, 0.3, 0.5 and 0.7 m.

Groundwater level was measured by a Dataqua DA-S-LRB 118 vented pressure transducer (Dataqua Elektronikai Kft., Balatonalmádi, Hungary) with time interval of 15 min. Manual groundwater level measurements were used to verify the reliability of the monitoring.

Water balance modeling with Hydrus 1-D

Model structure

The Hydrus 1-D ([34]) model was applied for the estimation of soil water content changes, actual transpiration and soil surface evaporation at both plots. Hydrus 1-D is a modeling software for analysis of water flow and solute transport in variably saturated porous media. The base of the model is the variable saturated vertical flow domain, where water flow is simulated.

The soil depth was set to three meters at both plots and was partitioned into seven soil layers according to the soil sampling (0-20 cm, 20-40 cm, 40-60 cm, 60-80 cm, 80-100 cm, 100-120 cm, 120-300 cm). The Hydrus model computed the variables of water flow for 200 soil horizons with gradually increments by depth. During the modeling daily time steps were used.

Daily potential transpiration T_P (mm day^-1) and evaporation E_P (mm day^-1) values together with the effective precipitation P_EFF (mm day^-1) constitute the upper boundary condition of the model. P_EFF is the difference between the gross precipitation (mm day^-1) and the interception loss (mm day^-1 - Fig. 3).

The net groundwater supply Q_NET (mm day^-1) defined the water flux at the lower boundary of the soil profile, which was computed using the high-frequency groundwater level data and the readily available specific yield S_y.

Water retention parameters and saturated hydraulic conductivity k_s (mm day^-1) were put into the Hydrus 1D program for each soil layers. The RETC software ([40]) was employed to get the soil water retention parameters from pairs of water content and pressure head data ([38]). Saturated hydraulic conductivity of the root zone was estimated for every soil layer using the method of Rawls et al. ([30]), although the final values were set by calibration. Root density was assigned to the profile based on field observations.

The actual values of transpiration T and soil surface evaporation E was computed using the potential values, soil hydraulic parameters, the available water content and the root density by the Hydrus 1-D model.

Effective precipitation

Effective precipitation was computed by estimation of the canopy and litter interception loss at the oak plot. The canopy interception loss was calculated applying the Gash model in the growing season (for details, see [9]). Besides the growing season, the daily interception method, proposed by Menzel ([23]) was used to estimate canopy interception loss. The method applies the daily sum of precipitation, the leaf area index and two species-specific parameters, describing the relation between LAI and the storage capacity ([43]) and the slope of the saturation curve ([29]).

Litter storage capacity was set during the model calibration.

Interception loss I was computed by a daily interception method at the fallow plot ([41] - eqn. 1):

where b is the soil cover fraction, which was set to unity since the soil surface was completely covered by vegetation. The empirical parameter a_i denotes the storage capacity of the vegetation by multiplying with the LAI (m² m^-2) and P is the daily sum of precipitation (mm day^-1).

Potential transpiration and soil surface evaporation

Potential transpiration T_P was defined as the canopy transpiration rate, when not limited by soil water deficits ([25] - eqn. 2):

where T_P is the potential transpiration (mm day^-1), K is the net shortwave radiation (MJ m^-2 day^-1), L is the net longwave radiation (MJ m^-2 day^-1), Δ is the slope of the saturation vapor pressure curve (kPa °C^-1), γ is the psychrometric constant (kPa °C^-1), VPD is the vapor pressure deficit (kPa), ρ_a is the air density (kg m^-3), c_a is the specific heat of moist air (kJ kg^-1 °C^-1), λ_v is the latent heat of vaporization (MJ kg^-1), ρ_w is the specific heat of air (kJ kg^-1 °C^-1), C_at is the aerodynamic conductance (m s^-1) and C_can is the canopy conductance (m s^-1). C_can was computed by multiplying the leaf area index by leaf conductance (m s^-1).

Leaf conductance C_leaf (m s^-1) was estimated after Jarvis ([15]) to express the effect of atmospheric demand on stomatal opening (eqn. 3):

where g_lmax and g_lmin are the maximum and minimum leaf conductance (m s^-1), respectively, and f_T, f_D, f_R, f_W, f_C are reduction factors, varying between zero and unity, that account for effects of temperature, vapor pressure deficit, radiation, leaf water stress and atmospheric carbon dioxide concentration on stomatal opening.

The reduction factors f_T and f_D were computed after Federer et al. ([7]) and Lohammar et al. ([20]). The radiation response f_R was estimated after Stewart ([36]) using the measured global radiation. The leaf water stress was omitted, thus f_W =1, and the effect of changing atmospheric carbon dioxide was also not considered, so f_C =1.

Potential soil surface evaporation E_P was also computed by the Penman-Monteith equation, based on the net radiation which reaches the soil surface ([31]).

Net groundwater supply

Diurnal fluctuation of groundwater levels occurs in shallow groundwater regions. The main inducing factor is the diurnal changes of transpiration during vegetation periods ([13]).

Gribovszki et al. ([12]) developed an empirical approach, based on the White-method ([44]) using the observed groundwater level data and the readily available specific yield S_y. The primary assumption was that Q_net displays a significant daily variation due to hydraulic gradient changes over the day. The maximum value of Q_net was calculated for each day by selecting the largest positive time-rate of change value in the groundwater level readings, while the minimum was obtained by computing the mean of the smallest time-rate of change in the predawn hours, by multiplying it by the value of S_y (eqn. 4):

where Q_NET is the net groundwater supply, h is the groundwater level, t is the time and S_y is the readily available specific yield.

The values of the minimum and maximum Q_net then were assigned to the locations of the groundwater level extremes. Spline interpolation was applied to interpolate the missing values between the minimum and maximum values. The daily sum of the net groundwater supplies was computed by summing up the half-hourly values.

Readily available specific yield was introduced by Meyboom ([24]) to denote the amount of water that can be released from the vadose zone during the time-frame of the daily fluctuations. S_y was estimated from soil texture data using the diagram of Loheide II et al. ([21]). The final values were set during the calibration of the model.

Actual transpiration and soil surface evaporation

Actual transpiration T was computed by integrating the actual transpiration rates of all model layers by the Hydrus 1-D model ([34]). The actual transpiration rate of a specific model soil layer is (eqn. 5):

where T_msl is the actual transpiration of one specific model soil layer, α(h) is the value of the S-shaped water uptake response function ([39]) and S_p is the potential uptake rate of that soil layer (mm day^-1), determined by the vertical root distribution.

Values of the parameter α(h) may be obtained as follows (eqn. 6):

where h is the pressure head (cm), h₅₀ represents the pressure head at which the water extraction rate is reduced by 50% (cm) and p is an experimental constant. The latter was set to three after van Genuchten ([39]) and h₅₀ was set to -4300 cm ([6]).

S_p is calculated using the value of daily potential transpiration and the normalized root water uptake distribution. For more details see Simunek et al. ([34]).

Actual soil surface evaporation E was estimated by the Hydrus model, based on soil water conditions and soil hydraulic properties using Darcy’s Law, bounded by the potential soil surface evaporation.

Model calibration

The observed soil water content and groundwater level data were employed to calibrate the Hydrus model at both plots. The model was calibrated specifically for this two year period and not used for future simulation.

During the model calibration, the storage capacity of vegetation and soil hydraulic parameters were modified to reduce the deviations between the modeled and observed soil water content of each layer and groundwater level.

At the oak plot, first the interception loss of the dormant season was calibrated. The storage capacity of stems and branches was set to 0.5 mm ([19]). Litter storage capacity was generally set to 0.5 mm. Storage capacity of the canopy in the vegetation period was calibrated to 1.17 mm at the maximum LAI 3.9 m² m^-2. A similar result (1.19 mm) was obtained by André et al. ([1]) in an oak dominated stand (LAI = 3.8 m² m^-2).

In the fallow plot, the best fit with the observed soil moisture was achieved by using a storage capacity of 1.0 mm during the vegetation period and 0.5 mm for the period of the dormant season.

During the calibration of the saturated hydraulic conductivity, the initially computed values were modified for the root-zone (<1.5 m at the oak and <0.8 m at the fallow plot) of the soil profile. The calibrated k_s values changed between 2.9 · 10³ mm day^-1 in the upper most horizon and 0.9 · 10³ mm day^-1 at about 1 m depth at the oak plot. In the saturated zone, the calibrated k_s value was around 1.4 · 10² mm day^-1 (Tab. 1).

At the fallow plot, the calibrated k_s values varied from 1.5 · 10² mm day^-1 to 1.1 · 10³ mm day^-1 in the root zone and around 0.8 · 10² mm day^-1 below the root zone.

Readily available specific yield was calibrated to 0.032 at the oak and to 0.029 at the fallow plot.

Model performance

The calibrated soil water content and groundwater levels were compared with measured values at both plots to check the performance of the model. The performance criteria were the coefficient of determination R², the root mean square error RMSE and the Nash-Sutcliffe modeling efficiency ME.

  Results and discussion 

Comparison of observed and modeled results

The calibrated soil water contents compared well with observations at both plots (Fig. 4).

The discrepancies at the beginning of the growing season of 2007 in measured vs. modeled soil water content may be attributed to the disturbance of the soil profile and vegetation at both plots during installation of the monitoring equipment. There were deviations at the fallow plot at a depth of 10 cm in December 2008 and January 2009, due to strong soil water freezing.

The calibrated groundwater levels compared quite well with measurements at both plots. Due to a malfunction of the pressure transducer, the continuous measurement failed at the fallow plot between November 2007 and June 2008. During this period regular manual groundwater depth measurements were employed to follow groundwater levels at the fallow plot (Fig. 5).

Tab. 3 presents the results of the performance criteria.

Comparison of the measured and calibrated soil water content and groundwater levels showed no systematic divergence at the plots (Fig. 6).

Water balance results

Total rainfall interception was twice as much in the forest than at the fallow plot considering the whole study period. In the 2007 growing season, 38% of the rainfall was intercepted at the oak plot, while at the fallow plot it was only 15% of the gross rainfall. As a consequence of more rainfall in 2008, the ratio of interception loss to precipitation decreased slightly at both plots.

According to the Hydrus model, the oak forest transpired approximately 33% more than the fallow vegetation, while groundwater consumption was three times higher during the study period.

Actual soil surface evaporation was only 4% of the total evapotranspiration at the oak plot, and 26% at the fallow plot during the whole study period. The low evaporation amount of the oak plot was the consequence of the high surface resistance due to the litter layer and the shading effect of the canopy. Available net radiation at the soil surface, and soil surface evaporation rates were higher at the fallow plot; thus the soil surface dried out more rapidly than at the oak plot.

Since weather conditions were quite contrasting during the growing seasons of 2007 and 2008, we decided to compare the water balances for both plots from 1^st of April until 30^th of September.

Groundwater consumption T_GW was computed using the temporal changes of the capillary fringe in the Hydrus model results, which was 66% of the total transpiration at the oak plot, while at the fallow plot it was only 38% in the growing season of 2007 (Fig. 7).

In this growing season, the water uptake showed the effect of drought in summer. Until the start of the summer, the transpiration from the unsaturated zone T_UZ was considerable from the total transpiration. The fraction of groundwater consumption increased to 60% at the fallow and to 90% at the oak plot in summer.

In the growing season of 2008 rainy weather provided a considerable amount of moisture for the unsaturated zone; thus the groundwater consumption was reduced from the preceding year. Groundwater consumption at the oak plot was 50% of the total transpiration, while at the fallow plot it was 25% in the growing season of 2008 (Fig. 8).

While in the growing season of 2007 both vegetation covers relied significantly on groundwater resources, in 2008 the evapotranspiration loss of groundwater was reduced considerably due to the rainy weather at both plots.

The water balance components of the growing seasons of 2007 and 2008 are shown in Tab. 4.

A meteorological tower was not set up at the oak plot, which contributed to the deviations of modeling results from the measurements. Air temperature and relative humidity at the fallow plot were applied above the canopy of the oak forest, which had an influence on the computation of potential transpiration and evaporation.

The empirical approach for computing Q_NET is further source of uncertainty, since it is sensitive to the quality, the time-step of the record of groundwater levels and the readily available specific yield. The groundwater level data was carefully examined and suspicious measurements were not considered for further computation. The time-step for calculation of Q_NET was half an hour, considering the 15 min frequency of groundwater level readings, suggested by Gribovszki et al. ([12]). The final value of the readily available specific yield was set to constant despite of its inherent variable nature in space and time.

  Conclusions 

Water balance components of an oak and fallow plot were estimated from April 1st 2007 to April 1st 2009 by calibrating the Hydrus 1-D model using soil moisture and groundwater level measurements. The study period included a dry (2007) and a wet growing season (2008).

For the entire study period, the Hydrus 1-D model results have shown that the total transpiration of the fallow plot was only two thirds of that obtained in the oak plot, while the soil surface evaporation in the oak plot was approximately one fifth of that in the fallow plot. The separation of transpiration into unsaturated transpiration and groundwater consumption has revealed that the groundwater consumption at the oak plot was almost three times higher than at the fallow plot. The groundwater consumption was close to 60% of the total transpiration at the oak forest and approximately 30% at the fallow plot.

By comparing the dry (2007) and wet (2008) growing seasons, we found that groundwater consumption was approximately 40% less in the wet than in the drier growing season, despite the fact that the groundwater level was deeper during the dry period. Thus, during the dry season both vegetation covers relied considerably on the available groundwater resources.

The results obtained of this investigation reinforce previous reports on higher groundwater consumption of forests as compared to other vegetation covers. Therefore, future afforestation in arid regions with shallow groundwater levels should pay attention to the large groundwater depleting effect of forests, especially in light of future climate change and human water extraction.

  Acknowledgments 

This research was funded by National R+D (NKFP 3B/2002/012 and NKFP 6-47/2005), EU Joint Development (TÁMOP-4.2.2-08/1-2008-0020 and TÁMOP 4.2.2.B-10/1-2010-0018 “Talentum”) projects and the HAS Bolyai scholarship. The authors would also like to acknowledge NASA for providing the MODIS data free of charge, and anonymous reviewers for providing input to improve the quality of the manuscript.

  References