One of the key inputs of a hydrologic budget is the potential evapotranspiration (PET), which represents the hypothetical upper limit to evapotranspirative water losses. However, different mathematical formulas proposed for defining PET often produce inconsistent results and challenge hydrological estimation. The objective of this study is to investigate the effects of the Priestley–Taylor (P–T), Hargreaves, and Penman–Monteith methods on daily streamflow simulation using the Soil and Water Assessment Tool (SWAT) for the southeastern United States. PET models are compared in terms of their sensitivity to the SWAT parameters and their ability to simulate daily streamflow over a five-year simulation period. The SWAT model forced by these three PET methods and by gauged climatic dataset showed more deficiency during low and peak flow estimates. Sensitive parameters vary in magnitudes with more skew and bias in saturated soil hydraulic conductivity and shallow aquifer properties. The results indicated that streamflow simulation using the P–T method performed well especially during extreme events’ simulation.

INTRODUCTION

Evapotranspiration (ET) constitutes a key element in developing strategies to optimize the use of water in a watershed and to introduce effective water management practices. It can be expressed as potential ET (PET; the amount of ET if ample water is available) or actual ET (AET). PET is a function of available energy, vapor pressure gradient, and vegetation type (Douglas et al. 2009) while AET is subject to spatial and temporal variability in soil type and moisture, water-table depth, and/or canopy characteristics. PET can be inferred through energy budget methods and/or hydrological water balance (WB) estimation (Samadi et al. 2014), with measurements of vegetation status and climate variables, i.e., stomatal conductance, temperature, solar radiation, air pressure, air humidity, and wind (Mu et al. 2011). PET is first simplified as a reference evapotranspiration (ET0), which is defined as the evaporation occurring from a land surface with ‘reference crop’ (e.g., green grass of 12 cm; Allen et al. 1998) on soil condition with sufficient water available.

Lu et al. (2005) reported that there are approximately 50 methods for PET estimation, including measurement-based, empirical, mass transfer-based, energy-based, temperature-based, and combination methods. Due to their input data requirements and different assumptions involved in the PET estimation methods, they provide different estimates. Direct measurements, such as those provided by lysimeters, pan evaporation method, eddy flux instrumentation or Bowen-ratio instrumentation, give point values, require constant attendance by skilled personnel, and are based on unverified assumptions (Morton 1983). Empirical PET techniques including temperature-based such as Hargreaves (HG) method, (Hargreaves et al. 1985), radiation-based such as Priestley–Taylor (P–T; Priestley & Taylor 1972) and combination of the two (e.g., Penman–Monteith (P–M); Monteith 1965) provide alternative measures of ET rates, and often offer promising results for the PET losses (e.g., Xu et al. 2006). Both radiation- and temperature-based PET models have been widely used for validating traditional hydrologic models that extrapolated findings at individual sites to the regional scale (e.g., Wilson et al. 2002). At present, empirical techniques have become more practicable in hydrological simulation since PET is a common input for watershed models that offers an upper limit to evapotranspirative water losses. Although a clear definition of the ‘best’ method for hydrologic computation is not still evident and the method choice is often subjective (e.g., Xu & Singh 2002; Douglas et al. 2009), primarily depending on the objectives of the study and the type of data available (Verstraeten et al. 2008).

A large body of literature related to the application of PET estimation methods exists in hydrology communities. For instance, Almorox & Grieser (2016) compared the Hargreaves–Samani (HS) method with the P–M equation (FAO-56) and calibrated the coefficients of HS for different climates as represented by the Köppen classification. Their results suggest that the correlation between long-term monthly means of HS and FAO-56 can be improved significantly by using climate-class specific coefficients. In another study, Koedyk & Kingston (2016) investigated the impact of error and uncertainty associated with PET on river flows using climate change scenario in the Waikaia River, New Zealand. They found less sensitivity of runoff projections to PET methods that needs to be explored more in the future. Alemayehu et al. (2015) studied the use of different sources of weather input data and PET estimation methods to derive ET in Kenya. They found that WB components such as ET, deep aquifer loss, and reevaporation from the shallow aquifer vary in magnitude, depending on the data and methods being used. Their study also showed that the HG method is more robust and realistic compare to the P–T and the P–M models.

Douglas et al. (2009) evaluated observed daily ET at 18 sites in Florida (USA) and stressed that the P–T performance is superior to the other methods (P–M and the Turc method (Tc)) for a variety of land covers. Weiß & Menzel (2008) compared the P–T method with the P–M and the HG methods, and advocated that the P–T results were closest to available pan evaporation data. Oudin et al. (2005) studied a total of 23 methods for PET estimation using a variety of micrometeorological input data. They compared the impact of the PET selection methods on four rainfall–runoff models over 308 watersheds located in France, Australia, and the United States. Their results showed that temperature-based PET estimates perform as well as or better than more physically based PET methods. Similar to Oudin et al. (2005), Vörösmarty et al. (1998) compared the performance of 11 different PET methods (simple to more physically based approaches) and achieved similar results. However, Oudin et al.'s (2005) study focused on systematic bias reduction using the P–M method prior to use in the rainfall–runoff models. They used relative performance of the PET methods in rainfall–runoff modeling and compared the computed PET with empirical ET measurements as subjective approaches.

This study focuses on the relative performance of the PET estimation methods in a distributed rainfall–runoff modeling and addresses the impact of the three different PET methods (i.e., P–T, HG, and P–M) on the watershed physical parameters over a forested wetland dominated ecosystem in the southeastern United States (SEUS). A few studies have used hydrology models to compare the relative performance of the PET selection methods in SEUS. For instance, Lu et al. (2005) compared several PET methods across 36 forested watersheds in SEUS with the result that the P–T method performed better than the rest of the approaches. Rao et al. (2011) explored different ways to identify appropriate PET models for two small forested watersheds in SEUS. They stressed that different PET methods provided significantly different estimates, but P–T resulted in the most reasonable estimates of forest ET rates.

While these studies provide insight into the PET estimation method, a challenge to selecting a representative PET method for heterogeneous regions is how the sensitivity of physical parameters may transfer/change depending on the method being used, and how (quantity) the sensitivity of physical parameters may affect the simulated flow regime. Due to scaling error in the PET estimation (e.g., Williams et al. 2004; Ford et al. 2007; Ngongondo et al. 2013) and often too simply derived as the residual of the WB (Domec et al. 2012), PET estimation methods have a large control on the sensitivity of physical parameters in rainfall–runoff modeling. However, although the impact of different PET estimation methods on the watershed physical parameters has been explored before (Andréassian et al. 2004; Zhao et al. 2013; Tamm et al. 2016; Singh et al. 2016; Koedyk & Kingston 2016), few studies have focused on the sensitivity of physical parameters to the PET estimation methods using uncertainty algorithm (e.g., Samadi et al. 2014; Alemayehu et al. 2015).

The objective of this study was to characterize the sensitivity of the Soil and Water Assessment Tool (SWAT) to the PET estimation methods using Bayesian uncertainty algorithm. The approach was to calibrate the three SWAT PET models (i.e., SWAT(P–T), SWAT(P–M), SWAT(HG)) using sequential uncertainty fitting (SUFI2) framework and assess the sensitivity of SWAT physical parameters during a simulation period over a heterogeneous watershed system in the SEUS. The results are compared in terms of their ability to capture the magnitudes of observed streamflow and flow duration curves (FDC). A unique aspect of this research is that this study is the first attempt to assess the sensitivity of physical parameters to the selected PET estimation methods in the SEUS.

The next section presents a description of the study area and introduces the SWAT model, the SUFI2 algorithm, and the PET methods. The methodology is then applied in the following section to identify sensitive parameters and to calibrate streamflow of a distributed rainfall–runoff model. The final section discusses and summarizes the findings of this research.

METHODOLOGY

Study area

The Waccamaw River watershed located in eastern North and South Carolina is a coastal plain watershed that was selected for this research (Figure 1). The study area is 3,116.85 sqkm (delineated by the SWAT model) and characterized by dense and rich vegetation with forested and non-forested wetland, wide and multiple floodplain, and low-gradient river in the downstream which comprises highly valued aquatic resources and is characterized by extensive alluvial riparian floodplains. In this landscape, hillslope processes dominate the hydrology of the upland, while hydrologic processes on relatively low-gradient poorly drained areas at downstream are usually dominated by shallow water table position. The river elevation gradients and sinuosity are very low, and do not vary significantly along its length. Wide floodplains occur along the entire length of the river, from headwaters to outlet.

Figure 1

The 311,685 ha Waccamaw River watershed used in SWAT. Meggett, Croatan, and Norfolk are the dominant soils on the floodplains. Overall, Rains and Meggett are the two major soil types in the Waccamaw River watershed.

Figure 1

The 311,685 ha Waccamaw River watershed used in SWAT. Meggett, Croatan, and Norfolk are the dominant soils on the floodplains. Overall, Rains and Meggett are the two major soil types in the Waccamaw River watershed.

This watershed typically experiences high climate variability and hydrological changes especially during wet and dry intervals. The climate of the watershed is humid subtropical in which winters are normally mild and summers are long and hot. Precipitation in the summer is dominated by convection storms, and in the winter by frontal boundaries. The long-term annual average temperature (1946–2009) is 16.88 °C, and averaged annual precipitation was a little more than 1,300 mm in the same period. Maximum temperature during the same period was 41.66 °C while minimum temperature was −14 °C.

The major soil type of this landscape is sandy loam with moderate permeability (see soil map in Figure 1). According to the SWAT output, the proportion of low storage shallow soil is around 90% which is restricted in moderately well drained (hydrologic group B soil) to poorly drained soil (hydrologic group D soil) with agriculture, rangeland, forested rangeland, and forested wetland land covers. Table 1 exhibits the percentage for each land use class (land use data of 2006) in the Waccamaw watershed. Forested wetland and evergreen forest are two major land use land cover classes in the study area.

Table 1

Land use area (sqkm) in the Waccamaw watershed. The area obtained using national land cover data of 2006

  Area
 
Land use type Land use description sqkm 
WATR Open water 48.71 1.6 
URLD Rural residential 112.81 3.7 
URMD Medium density residential 28.68 0.9 
URHD High density residential 4.40 0.1 
SWRN South western range 1.05 0.1 
FRSD Forest and woodland 4.21 0.1 
FRSE Evergreen forest 2.48 27.0 
FRST Mixed forest 843.97 1.4 
RNGB Range shrub land 42.20 10.8 
RNGE Grasslands/Herbaceous 334.82 5.0 
HAY Pasture 154.90 0.1 
AGRR Agricultural land-row crops 2.49 20.5 
WETF Forested wetland 641.75 27.7 
WETN Non-forested wetlands 864.11 1.0 
UIDU Industrial 28.99 0.1 
  Area
 
Land use type Land use description sqkm 
WATR Open water 48.71 1.6 
URLD Rural residential 112.81 3.7 
URMD Medium density residential 28.68 0.9 
URHD High density residential 4.40 0.1 
SWRN South western range 1.05 0.1 
FRSD Forest and woodland 4.21 0.1 
FRSE Evergreen forest 2.48 27.0 
FRST Mixed forest 843.97 1.4 
RNGB Range shrub land 42.20 10.8 
RNGE Grasslands/Herbaceous 334.82 5.0 
HAY Pasture 154.90 0.1 
AGRR Agricultural land-row crops 2.49 20.5 
WETF Forested wetland 641.75 27.7 
WETN Non-forested wetlands 864.11 1.0 
UIDU Industrial 28.99 0.1 

Sandy and sandy loam soils are the dominant types in the watershed. Approximately 67% of the soils are one of the four series, all of which are either hydrologic group B, D, or B/D. Only 7% of the soil is hydrologic group A; there is no group C soil (Table 2). Hydrologic group D soil (poorly drained) is adjacent to the main channel and hydrologic group B and B/D soils (mostly of the Rains series in the downstream reaches) are adjacent to those.

Table 2

Soil types in the Waccamaw watershed

 Hydrologic  Area
 
Series name soil group Drainage class sqkm 
Meggett Poorly drained 559.854 18.0 
Bohicket Very poorly drained 0.004 0.0 
Kenansville Well drained 49.424 1.6 
Leaf Poorly drained 144.951 4.7 
Woodington B/D Poorly drained 339.478 11.1 
Croatan Very poorly drained 473.114 15.2 
Rains B/D Poorly drained 614.282 19.7 
Norfolk Well drained 447.622 14.4 
Autryville Well drained 89.034 2.9 
Leon B/D Poorly drained 42.46 1.5 
Kureb Excessively well drained 15.614 0.5 
Trebloc Poorly drained 163.43 5.2 
Nansemond Moderately well drained 23.525 0.7 
Montevallo Well drained 35.07 1.1 
Osier Poorly drained 0.26 0.0 
Chisolm Well drained 118.71 3.8 
 Hydrologic  Area
 
Series name soil group Drainage class sqkm 
Meggett Poorly drained 559.854 18.0 
Bohicket Very poorly drained 0.004 0.0 
Kenansville Well drained 49.424 1.6 
Leaf Poorly drained 144.951 4.7 
Woodington B/D Poorly drained 339.478 11.1 
Croatan Very poorly drained 473.114 15.2 
Rains B/D Poorly drained 614.282 19.7 
Norfolk Well drained 447.622 14.4 
Autryville Well drained 89.034 2.9 
Leon B/D Poorly drained 42.46 1.5 
Kureb Excessively well drained 15.614 0.5 
Trebloc Poorly drained 163.43 5.2 
Nansemond Moderately well drained 23.525 0.7 
Montevallo Well drained 35.07 1.1 
Osier Poorly drained 0.26 0.0 
Chisolm Well drained 118.71 3.8 

Model calibration, validation, and sensitivity analysis was done using data from US Geological Service (USGS) gauging station at Longs. Research data were retrieved on a daily basis from the National Climatic Data Center and USGS portals on 6 September and 13 September 2012, respectively. Selected climate stations (Loris, Whiteville 7, and Longwood) are distributed evenly throughout the watershed.

ET processes in the coastal plain watershed

In the coastal plain ecosystem, surface condition is characterized by several features: the meteorological conditions in which the evaporation is taking place, degree of soil saturation, solar energy, vegetation and its seasonal dynamics, and the number of riparian buffers. Soil antecedent condition is different in the floodplain and hillslope in this ecosystem. More importantly, the surface is almost wet in riparian buffer whereas it is dry in uplands; therefore, evaporation/ET rates differ largely throughout the watershed system. Figure 2 visually illustrates the ET processes in the coastal plain watershed.

Figure 2

The effect of adverted air passing from dry condition over a wet condition in the coastal plain landscape.

Figure 2

The effect of adverted air passing from dry condition over a wet condition in the coastal plain landscape.

According to Monteith (1965), during ET processes water vapor is transported at a rate equal to the product of the water vapor content and the wind speed across a landscape. As shown in Figure 2, when air moves from a dried surface (upland) to a wet area (riparian buffer) in the coastal plain watershed, the concentration of water vapor increases at the transition to a higher value downwind. Dense vegetation and rich canopy as well as saturated/semi-saturated soil help to raise the concentration values throughout alluvial riparian floodplains. At the transition level, the AET rate may immediately maximize to a much higher value compared to that over the dry land, and then it may decline slowly to a value representative of the wetter condition. The lower AET over the dry surface in summer means the overpassing air may be hotter and drier, thus increasing the available heat energy to increase evaporation in the downwind wetter area (e.g., Morton 1983). These processes are a simple illustration of the AET processes in a coastal plain system and the context of this process is much more complicated in reality.

SWAT model

SWAT is a semi-distributed watershed hydrologic and water quality model initially developed by Arnold et al. (1993). The SWAT system is embedded within a geographic information system that can integrate various spatial environmental data including soil, land cover, climate, and topographic features. Based on this spatial information, a watershed is discretized into subwatersheds, and these in turn are subdivided into hydrological response units (HRUs). In a distributed hydrology model, a HRU is a unique area that is created by combination of land use, soil, and slope classes. These three classes have large effects on physical properties and runoff simulation in the model. Each process is treated at the HRU level on a t time step to provide aggregated inputs into a river network (Neitsch et al. 2011). SWAT then computes five different storages values at each HRU including snow, canopy storage, the soil profile with up to ten layers, a shallow aquifer, and a deep aquifer. SWAT is based on the WB equation (Equation (1)): 
formula
1
where ΔS is changes in water storage; t denotes daily time step; and P is the daily amounts of precipitation. Qtotal denotes the total amount of water yield. ET is evapotranspiration and Losses represents water losses via deep aquifer, percolation, and the evaporation from the shallow aquifer (REVAP).

SWAT offers three options for estimating PET using well-known empirical models, i.e., HG, P–T, and P–M. Users also can incorporate measured ET values in the model. A new CN method that includes a one-parameter soil moisture depletion curve to account for antecedent soil moisture conditions was used in this study to compute surface runoff. SWAT computes evaporation as a function of PET by combining ET values from soil and plants separately (e.g., Ritchie 1972). SWAT first evaporates any rainfall intercepted by the plant canopy. Next, it calculates the maximum amount of transpiration and sublimation/soil evaporation. The actual amount of sublimation and evaporation from the soil surface is then calculated. If snow is presented in the HRU, sublimation occurs and when there is no snow (such as this case study), only evaporation from the soil surface is calculated (Neitsch et al. 2005). The same procedure will be demonstrated for canopy and plant transpiration. Therefore, ET is the sum of the evaporation from the canopy and the soil as well as plant transpiration (Alemayehu et al. 2015).

PET estimation methods in the SWAT model

The P–M model accounts for the energy needed to sustain evaporation and for the mechanisms required to remove water vapor (e.g., Alemayehu et al. 2015). This method is usually used to estimate potential evaporation from a vegetated surface. PET can be estimated by the P–M method using the following formula (Equation (2)): 
formula
2
where PET represents potential evapotranspiration from a reference crop (mm/day); Δ is the slope of saturation vapor pressure–temperature curve (kPa/°C); denotes net radiation (MJ/m2/day); G is heat flux density to the ground (MJ/m2/day); is air density (kg/m3); denotes specific heat at constant pressure (MJ/kg/°C); represents saturation vapor pressure of air at height z (kPa); is water vapor pressure of air at height z (kPa); γ denotes psychrometric constant (kPa/°C); is plant canopy resistance or the (bulk) surface resistance (s/m); and is diffusion resistance of air layer (aerodynamic resistance; s/m).
The P–T model (Equation (3)) as a modification of the Penman model is a more theoretical equation, and can be used in areas of low moisture stress (Shuttleworth & Calder 1979). This method requires solar radiation, air temperature, and relative humidity. In practice, an empirical approximation of the Penman combination equation is made by the P–T to eliminate the need for input data other than radiation. It is reasoned that under ideal conditions, ET would eventually attain a rate of equilibrium for an air mass moving across a vegetation layer with an abundant supply of water, the air mass would become saturated and the actual rate of ET would be equal to the Penman PET rate. Under this circumstance, ET is referred to as equilibrium potential evapotranspiration (PETeq). The mass transfer term in the Penman combination equation approaches zero and the radiation terms dominate. Priestley & Taylor (1972) found that the AET from well-watered vegetation was generally higher than the equilibrium potential rate and could be estimated by multiplying the PETeq by a factor equal to 1.26: 
formula
3
where λ denotes the latent heat of vaporization (MJ/kg), Eo represents the potential evapotranspiration (mm d−1), is a coefficient, Δ indicates the slope of the saturation vapor pressure–temperature curve, de/dT (kPa/°C), γ signifies the psychrometric constant (kPa/°C), is the net radiation (MJ/m2/d), and G represents the heat flux density to the ground (MJ/m2/d). Although, the value of may vary throughout the day (Munro 1979), there is general agreement that a daily average value of 1.26 is applicable in a humid climate (Stewart & Rouse 1976; De Bruin & Keijman 1979; Shuttleworth & Calder 1979). Generally, the coefficient for an expansive saturated surface (like the coastal plain watershed) is usually greater than 1.0. This means that true equilibrium PET rarely occurs; there is always some component of advection energy that increases the actual evapotranspiration. Higher values of , ranging up to 1.74, have been recommended for estimating PET in more arid regions (ASCE 1990).
The HG method is empirical in nature and was developed based on daily grassed lysimeter data from Davis, California. Several improvements were made to the original equation (Hargreaves & Samani 1985) and the form used in SWAT was published in 1985 (Hargreaves et al. 1985). This method has a minimum weather data requirement, using only the maximum and minimum air temperature. The equation claims to self-compensate for the lack of solar radiation, humidity, and wind run data and to implicitly account for the effects of cloudiness, relative humidity, and vapor pressure deficit (Hargreaves & Allen 2003). The HG method equation is formulated as follows (Equation (4)): 
formula
4
where is incoming extraterrestrial solar radiation (MJ/m2/d). and are respectively, maximum and minimum air temperature for a given day (°C). When the HG method is used in SWAT, an adjusted evaporative water demand is first calculated by subtracting the free water evaporation from the canopy storage : 
formula
5
Next, the potential plant transpiration (Et) is estimated as a function of the simulated leaf area index (LAI): 
formula
6
The potential soil evaporation is a function of and the soil cover index : 
formula
7
The soil cover index depends on the above-ground biomass and residue (kg/ha). The actual soil evaporation is limited by the soil water content and is reduced exponentially when the soil water content drops below field capacity.

Model calibration and sensitivity analysis

This research examined a Bayesian framework (i.e., SUFI2) for the SWAT parameter sensitivity analysis. SUFI2 uses a multiple regression analysis to compute the statistics of parameter sensitivity. SUFI2 sensitivity analysis defines which parameters have significant influences on the model simulations by presenting the smallest P-value (values close to zero suggest high level of significance) and the largest t-stat (a measure of sensitivity, larger in absolute values are more sensitive). The t-stat represents the precision of the regression coefficient. It is the coefficient of a parameter divided by its standard error. If a coefficient is large compared to its standard error, then it is probably different from 0 and the parameter is sensitive. SUFI2 then compares the t-stat of a parameter with the values in the Student's t-distribution (defines the behavior of the mean of a sample) table to determine the P-value. The P-value for each term tests the null hypothesis that the coefficient is equal to zero (no effect). A low P-value (<0.05) indicates that the null hypothesis should be rejected. Thus, parameters with low P-values are likely to be meaningful additions to the model. Conversely, parameters with large P-values are more weakly associated with changes in the response (the parameter is not sensitive). Table 3 depicts the 15 parameters selected in the model calibration process for the SWAT model.

Table 3

Candidate parameters for SWAT calibration

SWAT parametera Name of SWAT parameter Absolute values Parameter type Identifier code 
CH_N2 Manning's ‘n’ value for the main channel −0.01–0.3 Channel route (.rte) Replaceb 
SOL_K Saturated hydraulic conductivity (mm/hr) 0–2,000 Soil (.sol) Relative 
CH_K2 Effective hydraulic conductivity in tributary channel alluvium (mm/hr) −0.01–500 Channel route (.rte) Replace 
SLSUBBSN Average slope length (m) 10–150 HRU (.hru) Replace 
ESCO Soil evaporation compensation factor 0–1 HRU (.hru) Replace 
OV_N Manning's ‘n’ value for overland flow 0.01–30 HRU (.hru) Replace 
SOL_BD Moist bulk density (Mg/m3 or g/cm30.9–2.5 Soil (.sol) Relative 
GWHT Initial groundwater height (m) 0–25 Groundwater (.gw) Relative 
EPCO Plant uptake compensation factor 0–1 HRU (.hru) Replace 
SHALLST Initial depth of water in the shallow aquifer (mm H2O) 1,000–50,000 Groundwater (.gw) Relative 
ALPHA_BF Base flow alpha factor (days) 0–1 Groundwater (.gw) Replace 
Curve number (CN) SCS runoff CN for moisture condition II 35–98 Groundwater (.gw) Relative 
RCHRG_DP Deep aquifer percolation fraction 0–1 Groundwater (.gw) Replace 
GW_REVAP Groundwater ‘revap’ coefficient 0.02–0.2 Groundwater (.gw) Replace 
GW_DELAY Groundwater delay time (days) 0–500 Groundwater (.gw) Replace 
SWAT parametera Name of SWAT parameter Absolute values Parameter type Identifier code 
CH_N2 Manning's ‘n’ value for the main channel −0.01–0.3 Channel route (.rte) Replaceb 
SOL_K Saturated hydraulic conductivity (mm/hr) 0–2,000 Soil (.sol) Relative 
CH_K2 Effective hydraulic conductivity in tributary channel alluvium (mm/hr) −0.01–500 Channel route (.rte) Replace 
SLSUBBSN Average slope length (m) 10–150 HRU (.hru) Replace 
ESCO Soil evaporation compensation factor 0–1 HRU (.hru) Replace 
OV_N Manning's ‘n’ value for overland flow 0.01–30 HRU (.hru) Replace 
SOL_BD Moist bulk density (Mg/m3 or g/cm30.9–2.5 Soil (.sol) Relative 
GWHT Initial groundwater height (m) 0–25 Groundwater (.gw) Relative 
EPCO Plant uptake compensation factor 0–1 HRU (.hru) Replace 
SHALLST Initial depth of water in the shallow aquifer (mm H2O) 1,000–50,000 Groundwater (.gw) Relative 
ALPHA_BF Base flow alpha factor (days) 0–1 Groundwater (.gw) Replace 
Curve number (CN) SCS runoff CN for moisture condition II 35–98 Groundwater (.gw) Relative 
RCHRG_DP Deep aquifer percolation fraction 0–1 Groundwater (.gw) Replace 
GW_REVAP Groundwater ‘revap’ coefficient 0.02–0.2 Groundwater (.gw) Replace 
GW_DELAY Groundwater delay time (days) 0–500 Groundwater (.gw) Replace 

aSWAT parameters are constructed based on SWACUP user manual (Abbaspour 2013).

bPrior distributions of aggregate parameters are based on Neitsch et al. (2011).

SUFI2 first defines the goal function with meaningful parameter ranges . Then, a Latin hypercube sampling is carried out in the initial hypercube to evaluate the corresponding goal function and to perform the sensitivity analysis (based on sensitivity matrix and the parameter covariance matrix). The 95PPU (calculated at the 2.5% and 97.5% levels of the cumulative distribution of output variables) of a sensitive parameter is then computed to find the best solutions (i.e., parameters which produce the optimal goal function). The 95PPU and model evaluation indices are then calculated. If model performances have satisfactory values, the parameter range travels from to as the posterior parameter distribution. Based on Bayes’ theorem, the probability density of the posterior parameter distribution will be driven from the prior density function to calculate the likelihood function of the model (see Yang et al. 2008). Uncertainty is then defined as discrepancy between measured and modeled variables where it is counted by variation between them.

Each SWAT model was calibrated using SUFI2 against five years of daily streamflow data (2003–2007), excluding the three-year warming-up period (i.e., 2000–2002). The calibration of each model was set for four iterations (each iteration representing 500 simulation runs). In each iteration, prior parameter ranges were updated by calculating the sensitivity matrix and the equivalent of a Hessian matrix, followed by the calculation matrix. Parameters were then updated in such a way that new ranges were always smaller than previous ranges and were centered on the best simulation.

Figure 3 illustrates a schematic framework for the SWAT calibration using the three PET estimation methods in SUFI2. SWAT first calculated the initial soil water content using effective rainfall, ET, surface runoff, the amount of water entering the vadose zone from the soil profile, and the amount of return flow. After saturating the vadose zone, PET is estimated and subtracted from runoff. Predicted runoff is routed to the channel to obtain the total streamflow (base flow and runoff). SWAT output was then linked with SUFI2 and examined sequentially (overall 2,000 runs) to calibrate streamflow. The best values of sensitive parameters were then used to optimize SWAT simulation and to obtain ET rates (as a part of WB component).

Figure 3

A conceptual representation of SWAT calibration for the Waccamaw watershed.

Figure 3

A conceptual representation of SWAT calibration for the Waccamaw watershed.

Model performance evaluation

In this study, the basic statistics used in the evaluation of the model performance are: Kling–Gupta efficiency (KGE; Gupta et al. 2009), Nash–Sutcliffe efficiency (NSE; Nash & Sutcliffe 1970), coefficient of determination R2 multiplied by the coefficient of the regression line (BR2), and sum of the squares of the differences of the measured and simulated values after ranking (SSQR).

KGE measure (Equation (8); 0 < =KGE < =1) is able to provide a diagnostically interesting decomposition of the NSE (Nash & Sutcliffe 1970) that allows to facilitate analysis of different elements (i.e., bias and correlation) in hydrologic simulation. Further, KGE keeps a high linear correlation between modeled and observed values (e.g., Gupta et al. 2009); hence the modeled value represents variability and better fit to high and low flow events. 
formula
8
where is the linear correlation coefficient between the observed flow and the simulated flow. is a measure of variability in the streamflow values (the standard deviation of the simulated flow over the standard deviation of the observed flow), and denotes the mean of simulated flow over the mean of observed flow (see Gupta et al. 2009).
The NSE (Equation (9)) is a normalized statistic that determines the relative magnitudes of the residual variance (‘noise’) compared to the measured data variance (Nash & Sutcliffe 1970). NSE ranges from −∞ to 1 where a high value (NSE > 0.5) indicates a good agreement between predicted and observed streamflow. 
formula
9
BR2 (0 < =bR2 < =1) is R2 (the coefficient of determination) multiplied by the coefficient of the regression line (Equation (10)). This function accounts for discrepancies in the magnitudes of two signals (depicted by b) as well as their dynamics (depicted by R2): 
formula
10
SSQR is similar to the mean square error method, and it is a measure to capture the fit of the frequency distribution to the observed and the simulated series. SSQR value of 0 indicates a perfect fit. The SSQ (Equation (11)) is defined as: 
formula
11

The SSQR is highly influenced by the magnitudes of the differences between the observed and simulated values. n is the number of pairs and i represents the rank. Overall, SSQR assesses the relative magnitudes of the bias in the model.

RESULTS AND DISCUSSION

Model application and sensitivity analysis

In this study, global sensitivity analysis was implemented for SWAT sensitivity test. This was conducted by allowing all parameters to be changed from their absolute values during the simulation period. This research evaluated and compared parameter sensitivity using a five-year calibration period with different starting points. Each simulation period was shifted by one year, such that subsequent periods have four years of data in common. Overall, ten different calibration periods were considered, and for each dataset, parameters were inferred using the SUFI2 model. This experience suggests that five years of daily streamflow data contain enough information about the parameters estimation of SWAT inputs using different PET estimation methods. Therefore, no significant variation in parameter sensitivity exists between different calibration periods for each PET setup model.

The sensitivity analysis revealed that when the PET estimation method changes, the sensitivity ranking of the parameters varies significantly. P–M method showed more sensitivity to ALPHA_BF, SLSUBBSN, and CH_K2 whereas HG is most sensitive to ESCO, CH_N2, and CN. Sensitivity results of the P–T method further demonstrated that the SWAT model is highly sensitive to SOL_K, CH_N2, and CH_K2 (see Table 4). Differences in sensitive parameters imply that the impact of the PET estimation methods on the ranking of the SWAT sensitive parameters is relatively high for the study watershed. In other words, changing the PET method has a great effect on physical parameters as well as streamflow simulation in this watershed. This could be attributed to the fact that the magnitudes and temporal dynamics of the streamflow showed a high sensitivity to physical parameters controlled by the PET values. Given the differences in the magnitudes of the estimated streamflow, SWAT parameters such as SOL_K, ESCO, and CN show a remarkable variation in rank by changing the PET estimation method. Because these parameters control physical processes and surficial fluxes in the watershed, these results stress that appropriate selection of the PET estimation method is important for appropriate representation of physical processes as well as estimating streamflow and WB components in the watershed under study.

Table 4

Parameter sensitivity results for the Waccamaw watershed using different PET estimation methods

 P–T
 
HG
 
P–M
 
SWAT parametera P-value t-stat Rank Optimal value P-value t-stat Rank Optimal value P-value t-stat Rank Optimal value 
v__CH_N2 7.48 0.10 4.19 0.101 0.144 −1.46 0.10 
r__SOL_K 0.07 1.83 0.058 0.31 0.31 0.51 −0.65 11 −0.88 
v__CH_K2 0.08 1.73 75.94 0.43 −0.78 93.67 0.1 1.62 185.32 
v__SLSUBBSN 0.17 1.37 71.32 0.138 1.48 74.31 2.65 77.54 
v__RCHRG_DP 0.21 −1.26 0.45 0.33 −0.96 0.60 0.17 1.35 0.49 
v__ESCO 0.25 −1.16 0.92 3.82 0.99 5.44 0.97 
v__OV_N 0.32 −1 0.56 0.41 −0.82 0.53 0.81 −0.23 14 0.26 
r__SOL_BD 0.33 0.97 0.002 0.94 0.063 15 0.001 0.15 1.43 0.60 
r__GWHT 0.41 0.83 1.33 0.46 0.73 10 0.66 0.96 −0.044 15 −0.50 
v__GW_REVAP 0.44 0.78 10 0.087 0.48 −0.7 11 0.18 0.142 −1.46 0.08 
v__EPCO 0.5 −0.63 11 0.55 −0.17 0.86 14 0.17 0.706 0.37 13 0.92 
r__SHALLST 0.67 0.43 12 −0.22 0.81 −0.22 13 0.97 0.47 0.72 10 −0.76 
v__ALPHA_BF −0.32 0.75 13 0.75 0.36 −0.9 0.6 0.13 −1.51 0.70 
r__CN −0.32 0.75 14 0.021 6.57 0.045 4.86 0.034 
v__GW_DELAY 0.94 0.07 15 207.71 0.7 0.37 12 307.99 0.69 0.4 12 249.86 
 P–T
 
HG
 
P–M
 
SWAT parametera P-value t-stat Rank Optimal value P-value t-stat Rank Optimal value P-value t-stat Rank Optimal value 
v__CH_N2 7.48 0.10 4.19 0.101 0.144 −1.46 0.10 
r__SOL_K 0.07 1.83 0.058 0.31 0.31 0.51 −0.65 11 −0.88 
v__CH_K2 0.08 1.73 75.94 0.43 −0.78 93.67 0.1 1.62 185.32 
v__SLSUBBSN 0.17 1.37 71.32 0.138 1.48 74.31 2.65 77.54 
v__RCHRG_DP 0.21 −1.26 0.45 0.33 −0.96 0.60 0.17 1.35 0.49 
v__ESCO 0.25 −1.16 0.92 3.82 0.99 5.44 0.97 
v__OV_N 0.32 −1 0.56 0.41 −0.82 0.53 0.81 −0.23 14 0.26 
r__SOL_BD 0.33 0.97 0.002 0.94 0.063 15 0.001 0.15 1.43 0.60 
r__GWHT 0.41 0.83 1.33 0.46 0.73 10 0.66 0.96 −0.044 15 −0.50 
v__GW_REVAP 0.44 0.78 10 0.087 0.48 −0.7 11 0.18 0.142 −1.46 0.08 
v__EPCO 0.5 −0.63 11 0.55 −0.17 0.86 14 0.17 0.706 0.37 13 0.92 
r__SHALLST 0.67 0.43 12 −0.22 0.81 −0.22 13 0.97 0.47 0.72 10 −0.76 
v__ALPHA_BF −0.32 0.75 13 0.75 0.36 −0.9 0.6 0.13 −1.51 0.70 
r__CN −0.32 0.75 14 0.021 6.57 0.045 4.86 0.034 
v__GW_DELAY 0.94 0.07 15 207.71 0.7 0.37 12 307.99 0.69 0.4 12 249.86 

aAggregate parameters are constructed based on SWACUP user manual (Abbaspour 2013).

‘v__’ and ‘r__’ means a replacement, and a relative change to the initial parameter values, respectively.

The fitted values for the SWAT parameters used in calibrations are shown in Figure 4. It is observed that ESCO is relatively lower for the P–M model, indicating that this method abstracts most of the evaporative demand from the lower soil layer. In contrast, both P–M and HG presented higher fitted values for ESCO, thus underestimated PET rates (higher value for the evaporation compensation factor also reveals decreasing in ET) compared to P–T. The fitted CN values are also higher for the HG and the P–M methods. SHALLST fitted values show that HG estimated this value very close to the SWAT upper absolute value while P–M underestimated and P–T provided median estimates. Regardless of hydraulic parameters, both soil and groundwater parameters show high variability in different PET methods. The fitted value for initial groundwater height is higher when the P–T method is employed, meaning that a large proportion of groundwater is retarded and stored in the soil. The lower pore space in a soil, the greater the value for the moist bulk density; this is a clear case for the HG fitted model with an overestimated SOL_BD. An increase in the soil moisture reduces the initial abstractions and, hence, higher runoff generation can be expected even at a lower CN value. As clearly observed from this result, the type of PET estimation method has great influences on the variability of physical parameters especially soil and groundwater properties.

Figure 4

Fitted calibrated SWAT parameters for different models: asterisks, optimized parameter values using SUFI2 algorithm; vertical columns, SWAT absolute parameter ranges.

Figure 4

Fitted calibrated SWAT parameters for different models: asterisks, optimized parameter values using SUFI2 algorithm; vertical columns, SWAT absolute parameter ranges.

Effect of PET estimation models on streamflow simulation

Figure 5 presents the measures of fit between the observed and SWAT simulated streamflow. Results revealed that the SWAT model with P–T as a reference PET method performed well compared to the rest of the methods. Based on Moriasi et al.’s (2007) qualitative rank, the SWAT(P–T) model showed very good performance while both SWAT(HG) and SWAT(P–M) showed satisfactory performances in simulating the observed daily streamflow during 2003–2007. This result is in agreement with the findings of Lu et al. (2005), Douglas et al. (2009), and Rao et al. (2011).

Figure 5

P–T model shows very good performance, while P–M and HG demonstrate satisfactory performances in simulating daily streamflow time series.

Figure 5

P–T model shows very good performance, while P–M and HG demonstrate satisfactory performances in simulating daily streamflow time series.

In addition, the error in streamflow prediction (SSQR; not shown in Figure 5) using different PET methods ranged from 61.43 (cubic meter per second; cm/s) to 271.55 cm/s and 294.57 cm/s, respectively, for P–T, HG, and P–M. HG and P–M methods yielded the highest error compared to the P–T method. This may be partly attributed to overpredicted/underpredicted ET during 2007 when the moisture was limited. The large discrepancy in 2007 was related to rainfall deficiency (>− 480 mm below average) and drought condition during the water year (1 October 2006 to 30 September 2007). However, spatial errors in rainfall data during the simulation period may also increase error and bias. For example, the average rainfall measured at the Loris station in 2006 was 1,253.8 mm compared to 1,580.3 mm at the Longwood gauge.

A close inspection of the simulated flow hydrographs reveals that the three SWAT models fairly capture the variability of the daily observed flows (Figure 6). The simulated streamflow follows the timing and magnitudes of the observed rainfall reasonably well, except during 2007 when this watershed experienced a severe drought condition. However, the SWAT(P–T) flow matches the rising limb and the recession curves of the observed hydrograph well compared to the rest of the models. A close observation also reveals that the SWAT(HG) model is better than the SWAT(P–M) model. This implies that when climate data are limited (radiation data are unavailable), the HG method is a proper replacement. This result is consistent with the findings of Lu et al. (2003, 2005) and Rao et al. (2011) for the coastal plain watershed.

Figure 6

Observed vs. simulated streamflow for the three SWAT PET methods.

Figure 6

Observed vs. simulated streamflow for the three SWAT PET methods.

Table 5 presents the quartiles statistics for the observed and simulated flows using the three PET methods. The simulated median flows using P–T, P–M, and HG are, respectively, 23.54 (cm/s), 27.82 (cm/s), and 28.42 (cm/s). Among them, the P–T model estimated median streamflow value relatively well. Based on this analysis, lower quantiles (low flow) were not successfully estimated by the three PET methods. In addition, high quantiles were mostly underpredicted by all methods, although Q95% and Q99% were approximately predicted by the P–T and the HG methods. Figure 7 presents the flow duration curves for the observed and simulated daily flows. The three SWAT models simulate slightly higher magnitudes during the dry season. All three methods also overestimated the lower band of FDC. The effect of the PET estimation methods seems to be significant during high and low flow estimates. In addition, the PET method has significant control on extreme events’ simulation.

Table 5

Streamflow statistics for different models

Quantiles Observed (cm/s) P–T (cm/s) HG (cm/s) P–M (cm/s) 
Q10% 1.11 10.83 16.08 10.04 
Q25% 8.01 14.92 21.09 18.72 
Q50% 23.13 23.58 28.42 27.82 
Q75% 46.16 47.46 43.65 44.65 
Q90% 85.40 84.32 79.78 75.26 
Q95% 124 110.58 105.06 96.51 
Q99% 173.70 163.90 175.98 149.86 
Quantiles Observed (cm/s) P–T (cm/s) HG (cm/s) P–M (cm/s) 
Q10% 1.11 10.83 16.08 10.04 
Q25% 8.01 14.92 21.09 18.72 
Q50% 23.13 23.58 28.42 27.82 
Q75% 46.16 47.46 43.65 44.65 
Q90% 85.40 84.32 79.78 75.26 
Q95% 124 110.58 105.06 96.51 
Q99% 173.70 163.90 175.98 149.86 
Figure 7

Observed and simulated flow duration curves with varying PET inputs; note that the vertical axis is log-transformed.

Figure 7

Observed and simulated flow duration curves with varying PET inputs; note that the vertical axis is log-transformed.

Effect of PET estimation method on AET estimation

Figure 8 depicts the SWAT-simulated AET, potential transpiration, potential soil evaporation, actual evaporation from the soil and the canopy interception, and the actual plant transpiration over the study area. Figure 8(a) shows monthly comparisons of the simulated AET using P–T, HG, and P–M methods. The differences in the simulated monthly AET while switching the PET estimation method could change AET amount to −2.2 mm/year (P–T to P–M), −6.7 mm/year (P–T to HG), and −4.45 mm/year (PM-HG), respectively. This indicates that the effect of the PET estimation method on the simulated AET is insignificant. On the other hand, the simulated AET values are less affected by the type of PET formula compared to PET. For all three formulas used in this study, higher simulated AET are obtained when the P–M method is selected whereas the HG and P–M methods presented moderate to lower AET values, correspondingly.

Figure 8

Monthly simulated AET based on the three PET methods (a). Annual potential and actual soil evaporation (b), potential and actual transpiration (c), as well as actual evaporation from canopy storage (d) based on P–T, P–M, and HG methods.

Figure 8

Monthly simulated AET based on the three PET methods (a). Annual potential and actual soil evaporation (b), potential and actual transpiration (c), as well as actual evaporation from canopy storage (d) based on P–T, P–M, and HG methods.

As expected, May to August is the peak time of ET in the Waccamaw watershed. It is interesting to note that AET was highly maximized during May to August 2005 with July having the highest AET rate (>277 mm; Figure 8(a)). Due to wet antecedent soil condition during 2003 (above normal rain of 1,580 mm) and somewhat during 2004 (normal rain of 1,240 mm but with wet antecedent condition), annual AET of 2005 was the lowest value in the study area. The annual AET remained close to the PET rate (around 90% of PET) for the years exceeding the long-term average rainfall and/or the years with just below the average but with the wet antecedent year. It appears in most of the year, the P–T method calculated moderate AET/PET values. In contrast, HG overestimated AET compared to P–T and P–M during summer season while underestimating during winter. The P–M method, on the other hand, overestimated AET during winter and underestimated during summer. From October 2006 to the end of 2007 with consistently below average monthly rainfall, all methods yielded high monthly AET with July 2007 being the highest value (>265 mm). A similar result was obtained by Lu et al. (2005), Douglas et al. (2009), and Rao et al. (2011). Wang et al. (2006) stressed that the sensitivity of the simulated AET to the choice of the PET estimation method were insignificant at a 5% significance level for the southeastern USA. Furthermore, Amatya & Jha (2011) advocated that SWAT is more sensitive to the choice of the PET estimation method; hence, it is important to select the most accurate and spatially stable model for the coastal plain watershed.

As shown in Figure 8(b) and 8(c), the magnitudes of the simulated potential plant transpiration and potential soil evaporation depend on the choice of the PET estimation method. Although potential plant transpiration estimation method is different in the PM, HG, and P–T methods, P–M estimated higher potential plant transpiration compared to the rest of the methods. Nevertheless, relatively lower is obtained from the simplified and less data-dependent method, i.e., the HG model. Overall, the annual average actual soil and plant transpiration reflects the effect of the choice of the PET estimation method on these physical properties, whereas the P–M method results in higher estimates in most cases (Figure 8(b)8(d)). On average, as shown in Figure 8(d), the P–T method demonstrated a higher annual evaporation from canopy storage, despite the relatively moderate PET estimates, as explained earlier. This is attributed to a higher EPCO value in the P–T method, as shown in Figure 4. The results for the EPCO parameter reveal the sensitivity of the simulated plant uptake compensation factor to the choice of the PET estimation method. Higher EPCO reveals that most of the transpiration demands are met by lower layers in the soil.

The violin plots of the marginal posterior distribution of the AET estimates are illustrated in Figure 9. A violin plot is a boxplot combined with kernel density plots, added on each side of the boxplot, to exhibit the probability distribution of the dataset. These plots were created using both PET and AET estimated values generated by the three PET estimation methods. PET estimates differ substantially in their posterior width. All models demonstrated similar median values for AET compared to PET. More dispersion and skewness is noted in the HG estimates. P–T provided longer probability tails (longer 5th and 95th percentile ranges) for both AET and PET estimates. The 25th and 75th percentile (thick black lines in Figure 9) ranges of PDF distribution are somewhat narrowed by P–M while the rest of models exhibited longer percentiles.

Figure 9

Violin plots of monthly ET and PET in each model. Thick black line and white dot show the 25th and 75th percentile ranges and median, respectively, and thin black line shows the 5th and 95th percentile ranges.

Figure 9

Violin plots of monthly ET and PET in each model. Thick black line and white dot show the 25th and 75th percentile ranges and median, respectively, and thin black line shows the 5th and 95th percentile ranges.

Spatial distribution of soil and groundwater physical properties

Since the selected PET method has a significant control on soil and groundwater parameters, this study analyzed the spatial distribution of those physical properties at subwatershed scale. Figures 1012 show the spatial distribution of soil and groundwater optimized contributions using P–T, HG, and P–M methods. The average annual soil saturation shows remarkable differences in the P–M method compared to the rest of the models. The P–M method overestimated both soil saturation and available water capacity compared to the rest of the methods. ET, on the other hand, was overestimated by the HG method, whereas wilting point was underestimated by the same method at each subwatershed scale. Both P–M and P–T estimated approximately similar wilting points because the structure of both methods is almost the same (P–T is a modification of P–M). Both methods are also more complex and more data dependent compared to HG. It seems wilting points are low at agricultural, hay and non-forested land covers, located mostly in subwatersheds 19, 24, and 26. Most agricultural plants will generally show signs of wilting lower and longer to restrict transpiration. However, soil types have significant influences on the wilting point value.

Figure 10

SWAT optimized values for soil and groundwater properties using the P–T method. P–T almost estimated similar values for wilting point, field capacity, and available water capacity at each subwatershed.

Figure 10

SWAT optimized values for soil and groundwater properties using the P–T method. P–T almost estimated similar values for wilting point, field capacity, and available water capacity at each subwatershed.

Figure 11

SWAT optimized values for soil and groundwater properties using the HG method. ET is overestimated whereas wilting point is underestimated.

Figure 11

SWAT optimized values for soil and groundwater properties using the HG method. ET is overestimated whereas wilting point is underestimated.

Figure 12

SWAT optimized values for soil and groundwater properties using the P–M method. Soil saturation and available water capacity were overestimated by this method.

Figure 12

SWAT optimized values for soil and groundwater properties using the P–M method. Soil saturation and available water capacity were overestimated by this method.

The field capacity is the amount of water remaining in the soil a few days after having been wetted and after free drainage has ceased. It appears this value is overestimated by the P–M method while HG underestimated it compared to P–T. The total available water (holding) capacity is the proportion of water available, stored, and/or released between field capacity and the wilting point water contents. The soil types with higher total available water content are generally more conducive to high biomass productivity because they can provide adequate moisture to plants during water shortage. Soils with low available water capacity may need more irrigation and this is the case for shallow rooted crops (i.e., potatoes) located on sandy soil in the study area. Higher available water capacity leads to increased soil moisture in the soil profile, which in turn reduces the initial abstractions, hence, higher surface runoff generation.

In addition, analysis revealed that the highest ET rate is calculated for subwatersheds 25 and 5 with forested wetland land cover due to its potential for ET processes while the least ET is calibrated, respectively, for range and agricultural lands at subwatersheds 26 and 19. Spatial distribution of ET estimation methods also indicates that the HG model overestimated ET for forested wetland, agricultural, and rangeland land covers compared to the rest of the methods (for instance, see HG estimated values for subwatersheds 20 and 25). All three methods also calculated the same values for the forest land cover. Overall, ET distribution showed that agricultural land covers presented moderately higher ET values. Furthermore, all methods presented approximately similar values for groundwater contribution at each subwatershed level.

CONCLUSIONS

The aim of this study was to examine the sensitivity of SWAT physical parameters to changes in the PET methods over a coastal plain ecosystem in the southeastern USA. The output of SWAT was compared in terms of parameter sensitivity and streamflow simulations using Bayesian uncertainty algorithm.

The use of ET equations in the SWAT model cause different degrees of sensitivity (and consequently increasing bias) for physical parameters because of either the assumptions made in interrelationship between climatic variables associated with the PET equation or simplified relationship between atmospheric variables, soil and plant/crop computations. Analysis shows that among all three PET estimation models examined here, radiation-based P–T equation which was developed for warm and humid climate condition, performed best for the Waccamaw coastal plain watershed when judged by streamflow simulation, sensitivity analysis, and physical parameter estimation. This reveals that the radiation-based method is efficient (solar radiation, air temperature, and relative humidity control ET process); hence the P-T model may be preferable to the HG and P–M models for the coastal plain with wet and humid surface, as stressed by Lu et al. (2005), Wang et al. (2006), and Rao et al. (2011).

Although there are large differences in PET estimates by various methods during the dormant season (winter) compared to the growing season (spring and summer), modeling results revealed insignificant differences among these three methods in respect to AET estimates. Further, PET estimates remarkably differed during peak and low flow events as well as winter and summer seasons. This analysis reveals that low and high flows are more sensitive to the changes in the PET estimation method. Thus, a careful selection of the PET method is required when using SWAT for extreme simulation (drought and flood) in the study watershed.

The adequacy of the assumptions made in the P–T equation for the Waccamaw watershed has been validated by a number of independent investigations before (Lu et al. 2005; Wang et al. 2006; Samadi et al. 2014) in which it was commonly found that, in vegetated areas with no or very small water deficiency, approximately 95% of the annual evaporative demand is supplied by radiation (Stagnitti et al. 1989). This research finding is also consistent with the study by Rao et al. (2011), who compared three common PET models (FAO-56 grass reference ET, Hamon PET, and P–T) in the southeastern Appalachian mountain region. They found that the radiation-based equation (i.e., P–T method) generally performed better than those that included only temperature-related input variables.

Overall, this results of this study showed that AET represents 35%–53% of precipitation losses in hydrologic modeling of the Waccamaw coastal plain watershed. It is fairly in agreement with Palmroth et al. (2010), who stressed that 30% of AET was explained by the variation in precipitation in the coastal plain watershed. However, the quality of input data and methods for calculating PET (and hence, AET) matters and should receive more attention when applying hydrology models in data-scarce coastal plain watersheds. It is recommended that a comparison of SWAT and remotely sensed ET estimates using climate change scenarios (e.g., Awan et al. 2016; Tamm et al. 2016) be made to provide spatial distributed ET values at subwatershed or HRU level.

ACKNOWLEDGEMENTS

The author appreciates those people and agencies that assisted in accessing research data and information. This research was partially supported by South Carolina Sea Grant Consortium (Grant 15520-GA11) and the University of South Carolina (Grant 15520-16-40787). I also wish to thank the anonymous reviewers and the editor for their valuable comments.

REFERENCES

REFERENCES
Abbaspour
K. C.
2013
SWAT-CUP 2012: SWAT Calibration and Uncertainty Programs–A User Manual
.
Swiss Federal Institute of Aquatic Science and Technology
,
Eawag, Duebendorf
,
Switzerland
, p.
103
.
Allen
R. G.
Periera
L. S.
Raes
D.
Smith
M.
1998
Crop Evapotranspiration: Guidelines for Computing Crop Requirements, Irrigation and Drainage Paper No. 56. FAO, Rome, p. 300
.
American Society of Civil Engineers, Committee on Irrigation Water Requirements of the Irrigation and Drainage Division of the ASCE
1990
Evapotranspiration and Irrigation Water Requirements: a Manual
, p.
332
.
Arnold
J. G.
Allen
P. M.
Bernhardt
G.
1993
A comprehensive surface-groundwater flow model
.
J. Hydrol.
142
,
47
69
.
Awan
U. K.
Liaqat
U. W.
Choi
M.
Ismaeel
A.
2016
A SWAT modeling approach to assess the impact of climate change on consumptive water use in Lower Chenab Canal area of Indus basin
.
Hydrol. Res
.
47
(
5
),
1025
1037
.
DOI: 10.2166/nh.2016.102
.
De Bruin
H. A. R.
Keijman
J. Q.
1979
The Priestley-Taylor evaporation model applied to a large, shallow lake in the Netherlands
.
J. Appl. Meteor.
18
,
898
903
.
Douglas
E. M.
Jacobs
J. M.
Sumner
D. M.
Ray
R. L.
2009
A comparison of models for estimating potential evapotranspiration for Florida land cover types
.
J. Hydrol.
373
(
3–4
),
366
376
.
Ford
C. R.
Hubbard
R. M.
Kloeppel
B. D.
Vose
J.
2007
A comparison of sap flux-based evapotranspiration estimates with catchment-scale water balance
.
Agr. Forest Meteor.
145
,
176
185
.
Hargreaves
G. H.
Allen
R. G.
2003
History and evaluation of Hargreaves evapotranspiration equation
.
J. Irrig. Drain. Eng.
129
(
1
),
53
63
.
Hargreaves
G. H.
Samani
Z. A.
1985
Reference crop evapotranspiration from temperature
.
Appl. Eng. Agric.
1
(
2
),
96
99
.
Hargreaves
G. L.
Hargreaves
G. H.
Riley
J. P.
1985
Agricultural benefits for Senegal River watershed
.
J. Irrig. Drain. Eng.
111
(
2
),
113
124
.
Koedyk
L. P.
Kingston
D. G.
2016
Potential evapotranspiration method influence on climate change impacts on river flow: a mid-latitude case study
.
Hydrol. Res
.
47
(
5
),
951
963
.
DOI: 10.2166/nh.2016.152
.
Lu
J.
Sun
G.
McNulty
S. G.
Amatya
D. M.
2003
Modeling actual evapotranspiration from forested watersheds across the southeastern United States
.
J. Am. Water Resour. Ass.
39
,
886
896
.
Lu
J.
Sun
G.
McNulty
S. G.
Amatya
D. M.
2005
Comparison of six potential evapotranspiration methods for regional use in the southeastern United States
.
J. Am. Water Resour. Ass.
41
,
621
633
.
Monteith
J. L.
1965
Evaporation and environment
. In:
The State and Movement of Water in Living Organisms, Symposium Society Experimental Biology
(
Fogg
G. E.
, ed.).
Cambridge University Press
,
London
,
UK
, p.
432
.
Moriasi
D. N.
Arnold
J. G.
Liew
M. W.
Van Bingner
R. L.
Harmel
R. D.
Veith
T. L.
2007
Model evaluation guidelines for systematic quantification of accuracy in watershed simulations
.
Trans. ASABE
50
(
3
),
885
900
.
Mu
Q.
Zhao
M.
Running
S.
2011
Running improvements to a MODIS global terrestrial evapotranspiration algorithm
.
Remote Sens. Environ.
115
,
1781
1800
.
Munro
D. S.
1979
Daytime energy exchange and evaporation from a wooded swamp
.
Water Resour. Res.
15
(
5
),
1259
1265
.
Nash
J. E.
Sutcliffe
J. V.
1970
River flow forecasting through conceptual models: Part 1. A discussion of principles
.
J. Hydrol.
10
(
3
),
282
290
.
Neitsch
S. L.
Arnold
J. G.
Kiniry
J. R.
Williams
J. R.
2005
Soil and Water Assessment Tool Theoretical Documentation Version 2005
.
Grassland Soil and Water Research Laboratory, Agricultural Research Service
,
Texas
,
USA
, p.
618
.
Neitsch
S. L.
Arnold
J. G.
Kiniry
J. R.
Williams
J. R.
2011
Soil and water assessment tool theoretical documentation version 2009
.
Tech. Rep. 406
,
Texas Water Resources Institute, Texas A&M University, College Station
,
Texas
,
USA
, p.
647
.
Palmroth
S.
Katul
G. G.
Dafeng
H.
McCarthy
H. R.
Jackson
R. B.
Oren
R.
2010
Estimation of long-term basin scale evapotranspiration from streamflow time series
.
Water Resourc. Res.
46
,
W10512
.
Priestley
C. H. B.
Taylor
R. J.
1972
On the assessment of surface heat flux and evaporation using large-scale parameters
.
Mon. Weather Rev.
100
(
2
),
81
92
.
Rao
L. Y.
Sun
G.
Ford
C. R.
Vose
J. M.
2011
Modeling potential evapotranspiration of two forested watersheds in the southern Appalachians
.
Trans. ASABE
54
(
6
),
2067
2078
.
Ritchie
J. T.
1972
Model for predicting evaporation from a row crop with incomplete cover
.
Water Resour. Res.
8
(
5
),
1204
1213
.
Samadi
S.
Tufford
D.
Carbone
G.
2014
Improving hydrologic predictions of distributed watershed model via uncertainty quantification of evapotranspiration methods
. In:
South Carolina Water Resource Conference
,
15–16 October 2014
,
Columbia, SC
,
USA
.
Shuttleworth
W. J.
Calder
I. R.
1979
Has the Priestley-Taylor equation any relevance to forest evaporation
?
J. Appl. Meteor.
18
(
5
),
639
646
.
Singh
H. V.
Kalin
L.
Morrison
A.
Srivastava
P.
Lockaby
G.
Pan
S.
2016
Post-validation of SWAT model in a coastal watershed for predicting land use/cover change impacts
.
Hydrol. Res.
46
(
6
),
837
853
.
Stagnitti
F.
Parlange
J. Y.
Rose
C. W.
1989
Hydrology of a small wet catchment
.
Hydro. Proc.
3
,
137
150
.
Stewart
R. B.
Rouse
W. R.
1976
A simple method for determining the evaporation from shallow lakes and ponds
.
Water Resour. Res.
12
(
4
),
623
628
.
Tamm
O.
Luhamaa
A.
Tamm
T.
2016
Modelling future changes in the North-Estonian hydropower production by using SWAT
.
Hydrol. Res
.
47
(
4
),
835
846
.
DOI: 10.2166/nh.2015.018.
Verstraeten
W. W.
Veroustraete
F.
Feyen
J.
2008
Assessment of evapotranspiration and soil moisture content across different scales of observation
.
Sensors
8
,
70
117
.
Wang
X.
Melesse
A. M.
Yang
W.
2006
Influences of potential evapotranspiration estimation methods on SWAT's hydrologic simulation in a northwestern Minnesota watershed
.
Am. Soc. Agric. Biol. Engr.
49
(
6
),
1755
1771
.
Williams
D. G.
Able
W. C.
Ultine
K. H.
Hoedjes
J. C. B.
Yepez
E. A.
Immoneaux
V. S.
Er-raki
S.
Oulet
G. B.
de Bruin
H. A. R.
Chehbouni
A. C.
Hartogensis
O. K.
Timouk
F.
2004
Evapotranspiration components determined by stable isotope, sap flow and eddy covariance techniques
.
Agr. Forest Meteor.
125
,
241
258
.
Wilson
K.
Goldstein
A.
Falge
E.
Aubinet
M.
Baldocchi
D.
Berbigier
P.
Bernhofer
C.
Ceulemans
R.
Dolman
H.
Field
C.
Grelle
A.
Ibrom
A.
Law
B. E.
Kowalski
A.
Meyers
T.
Moncrieff
J.
Monson
R.
Oechel
W.
Tenhunen
J.
Verma
S.
Valentini
R.
2002
Energy balance closure at FLUXNET sites
.
Agr. Forest Meteor.
113
,
223
243
.
Xu
C.-Y.
Singh
V. P.
2002
Cross-comparison of mass-transfer, radiation and temperature based evaporation models
.
Water Resour. Manage.
16
,
197
219
.
Yang
J.
Abbaspour
K. C.
Reichert
P.
Yang
H.
2008
Comparing uncertainty analysis techniques for a SWAT application to Chaohe Basin in China
.
J. Hydrol.
358
,
1
23
.
Zhao
L. L.
Xia
J.
Xu
C.-Y.
Wang
Z. G.
Sobkowiak
L.
Long
C. R.
2013
Evapotranspiration estimation methods in hydrological models
.
J. Geogr. Sci.
23
(
2
),
359
369
.