Using the SWAT model to improve process descriptions and define hydrologic partitioning in South Korea

Watershed-scale modeling can be a valuable tool to aid in quantification of water quality and yield; however, several challenges remain. In many watersheds, it is difficult to adequately quantify hydrologic partitioning. Data scarcity is prevalent, accuracy of spatially distributed meteorology is difficult to quantify, forest encroachment and land use issues are common, and surface water and groundwater abstractions substantially modify watershed-based processes. Our objective is to assess the capability of the Soil and Water Assessment Tool (SWAT) model to capture eventbased and long-term monsoonal rainfall–runoff processes in complex mountainous terrain. To accomplish this, we developed a unique quality-control, gap-filling algorithm for interpolation of high-frequency meteorological data. We used a novel multi-location, multi-optimization calibration technique to improve estimations of catchment-wide hydrologic partitioning. The interdisciplinary model was calibrated to a unique combination of statistical, hydrologic, and plant growth metrics. Our results indicate scale-dependent sensitivity of hydrologic partitioning and substantial influence of engineered features. The addition of hydrologic and plant growth objective functions identified the importance of culverts in catchment-wide flow distribution. While this study shows the challenges of applying the SWAT model to complex terrain and extreme environments; by incorporating anthropogenic features into modeling scenarios, we can enhance our understanding of the hydroecological impact.


Introduction
Land use and land cover (LULC) distribution can have a substantial influence on catchment water balance due to localized precipitation, evaporation, transpiration, soil moisture redistribution, and crop associated temporal variations in surface runoff.The effects of land use change, including deforestation (Forti et al., 1995), agricultural intensification (Berka et al., 2001), yearly variations in agricultural land use (Tilman et al., 2002), and construction of roads, culverts, and sediment detention ponds (Strauch et al., 2014) on stream discharge and water quality occur at many spatial and temporal scales.Deforestation significantly affects streamflow Published by Copernicus Publications on behalf of the European Geosciences Union.
characteristics (Calder, 1992) by increasing erosion and decreasing soil moisture and soil nutrient concentrations.Agricultural intensification influences surface runoff by altering infiltration, evaporation, and timing of runoff.As agricultural land use increases, the need for water resources management increases, particularly in complex topography driven by extreme events.
The water resources of the Haean catchment in South Korea are important to quantify because the catchment represents an important contributor to the Han River and the Soyang Lake watershed, which is a major drinking water source for major metropolitan areas including the city of Seoul (Jo and Park, 2010).The catchment is also a significant source of sediment and nutrients due to the high agricultural activity and forest encroachment (Jung et al., 2012;Lee et al., 2014).Small-scale agriculture is the largest economic activity within the basin, engaging 85 % of the population and up to 44 % of the available land area within the catchment.Increasing agricultural encroachment into the forest region imposes a significant risk to water yield and quality with a reduction in forested area by 37 % over the past 20 yr (Kim et al., 2011).Furthermore, routing and flow management in Haean has significantly increased the erosive power and decreased infiltration during individual events (Arnhold et al., 2013).Previous studies have suggested an appreciable decline in aquatic species, attributed in large part to an increase in fine grain sediment erosion and nutrient concentrations (B.Kim, personal observation, 2010;Jun, 2009).Since the end of the Korean War in 1953, a variety of amelioration measures such as river regulation, installation of catchment drainage systems, and waste water treatment plants (WWTPs) have been implemented in order to enlarge communities and increase local agricultural production.These measures have led to a change in the catchment-wide water balance, spatiotemporal nutrient dynamics, and floodplain ecology (Jun, 2009).Several conservation projects have been implemented within the Haean catchment and throughout South Korea to limit and effectively manage soil erosion including retention pond construction, modification of riparian channel widths, and channel reinforcement.Consequently, the landscape has been intensively altered, creating a mosaic of ecohydrologic landscape patterns.Surface water and groundwater abstractions, dam and reservoir operations, and engineered hydraulic structures (culverts, sediment ponds, and roads) have disrupted the natural hydrology of the catchment.In higher elevations, surface water flow has been observed to be entirely depleted over extended stretches due to domestic and irrigation abstractions for dryland farms (Shope et al., 2013).Previous research has indicated that seasonal precipitation, as well as individual events, influences the hydrologic flushing of organic materials from the land surface (Jung et al., 2012;Lee et al., 2014).The longterm interdisciplinary research group TERRECO (Tenhunen et al., 2011), has collected spatiotemporal terrestrial surface runoff measurements to calculate sediment yield (Arnhold et al., 2013), conduct dye tracer experiments to estimate soil structure and variably saturated flow and transport processes (Ruidisch et al., 2013), and examine groundwater and surface water exchange on spatiotemporal fluxes and nearstream biogeochemistry (Bartsch et al., 2014).To quantify overland runoff, sediment transport, and soil loss from individual crops under specific management practices, it is critical to understand sustainable resource allocation and scenario implications in this agriculturally productive, complex terrain.
Coupled hydrological and crop production watershedscale models are a useful tool to simulate the interactions of catchment physical characteristics, agricultural practices, and weather inputs on the water yield and to evaluate conservation practices in locations with limited observational data (Cho et al., 2012).Model scenarios can be helpful in identifying reasonable measures for assessing environmental ecological status (Lam et al., 2012;Volk et al., 2009).Gassman et al. (2007) found that the distributed Soil and Water Assessment Tool (SWAT) model was a promising model for predominately agricultural watersheds located throughout the world when compared to several other integrated watershed models.SWAT has also been successfully applied in a wide variety of data-limited studies, particularly in South Korea (Lee et al., 2012(Lee et al., , 2011;;Stehr et al., 2008;Mekonnen et al., 2009).We use the SWAT model because it is a welldocumented, efficient model that couples long-term climate, land use, and management practices to evaluate catchmentwide hydrology.
This study builds upon multiple research investigations distributed throughout the Soyang Lake watershed by implementing the SWAT ecohydrologic model within the Haean catchment to quantify hydrologic processes and catchmentwide flow partitioning.Our objectives are to (1) assess the potential of a spatiotemporal algorithm to improve discretization of monitored precipitation, (2) characterize the spatiotemporal river discharge patterns at multiple locations throughout the monsoon driven catchment through multiobjective optimization, (3) determine the capability of the SWAT model to capture daily monsoonal rainfall-runoff processes in complex mountainous terrain, and (4) quantify the significance of engineered structures (roads, culverts, sedimentation ponds) on flow partitioning.To accomplish these objectives, we utilized robust and comprehensive, spatiotemporal river discharge estimates at 14 locations throughout the Haean catchment to quantify flow partitioning.We discuss the construction of the ecohydrologic SWAT model for the Haean catchment, the selection and sensitivity of model parameters, and the calibration and validation of the model.Finally, we evaluate three different river routing systems including (1) the surface water drainages; (2) a combination of the rivers and engineered culverts; and (3) the rivers, culverts, and road network, to identify flow partitioning throughout the catchment.

45
Figure 1-Haean study area within the Lake Soyang watershed is located in northeastern South Korea along the demilitarized zone (DMZ) border with North Korea.The regional KMA weather station and local meteorological stations are denoted with white circles and (WS).
River discharge monitoring locations are denoted by (S) and the yellow squares.

Catchment characteristics
The Haean catchment study area (38.239-38.329• N,  128.083-128.173• E) is located in the Gangwon Province of the northeastern portion of South Korea along the demilitarized zone (DMZ) between South and North Korea (Fig. 1).The 62.7 km 2 catchment has a unique bowl-shaped physiographic characteristic with elevation ranging between 339 to 1321 m a.s.l., which drastically alters the local meteorological conditions.The catchment drainage is the Mandae River with a maximum length of 8.6 km.Limited historical observations are available, although this is typical for most areas outside of Europe and North America.The average catchment discharge at the outlet is 4.32 m 3 s −1 (1.20-379 m 3 s −1 ) while the average discharge at the S1 headwater monitoring location is 0.03 m 3 s −1 (1.4 × 10 −4 -10.0 m 3 s −1 ).The catchment hydrology is further described in Shope et al. (2013).The catchment is 56 % forested and 44 % agricultural LULC.
Geologically, the basin is composed of a Precambrian gneiss complex at the higher elevation mountain ridges and a highly weathered Jurassic biotite granite intrusion that was subsequently eroded throughout the central portion of the catchment (Kwon et al., 1990).Alluvium generally extends up to 2 m in depth and bedrock is typically observed between 20 and 45 m below land surface in the catchment interior.Surficial soil texture is typically saprolitic sand and sandy loam with high infiltration capacity (Arnhold et al., 2013;Jo and Park, 2010).
The climate in South Korea is humid continental to humid subtropical, influenced by the East Asian summer monsoon and early autumn typhoons.The monsoon season extends from the end of June through the end of July, followed by scattered events through early September, with up to 70 % of the total annual precipitation between the months of June and August.The average annual rainfall over the most recent 12 yr of record is 1514 mm (930 to 2299 mm yr −1 ) with a maximum precipitation as high as 48.6 mm h −1 or up to 223.2 mm d −1 .The average annual temperature is 8.65 ± 0.35 • C ranging between −26.9 • C in January to 33.4 • C in August.Choi et al. (2010) found that the temperature lapse rate within the Haean catchment ranged between −0.56 • C 100 m −1 throughout the spring to +1 • C 100 m −1 during early morning inversions after many consecutive sunny days.

Model description
The SWAT model is a continuous, physically based, distributed model originally developed to predict the long-term impact of climate and land use management practices on hydrologic, sediment, and agricultural chemical yields in large, complex basins (Arnold et al., 1998).Essentially, SWAT uses the water balance approach to simulate watershed hydrologic partitioning as described by Neitsch et al. (2010).Catchments are divided, typically on a topographic basis, into spatially linked subbasins and the subbasins are segregated into unique hydrological response units (HRUs) by integrating the combination of LULC, soil type, and slope to describe the system physical heterogeneity.The modeled hydrological components include surface runoff, percolation, lateral flow, groundwater flow, evapotranspiration (ET), and transmission losses.The simulation of watershed hydrology with SWAT is split into the land phase and the channel or routing phase of the hydrologic cycle, which controls the amount of water, sediment, and nutrients into the main channel in each C. L. Shope et al.: Landscape complexity and ecosystem modeling with the SWAT model subbasin and through the channel network to the watershed outlet (Neitsch et al., 2011).Incoming precipitation is partitioned into canopy storage, infiltration, and surface runoff through either the SCS (Soil Conservation Service) curve number (CN) method (U.S. D.A., 1972) or the Green-Ampt (Green and Ampt, 1912) method.Daily runoff volume from the SCS retention parameter can be calculated through the shallow soil water content or through accumulated plant ET.The SCS curve number method with calculated plant evapotranspiration was selected for the Haean catchment simulations.The hydrologic condition of the vegetation is important in determining CN for individual HRUs (U.S.D. A., 1972).Therefore, the distributed CN was further modified within individual HRUs through time-variable LULC characterization and crop growth.The model uses the modified Rational Method to estimate peak flow (Neitsch et al., 2011).Runoff in SWAT is aggregated from the HRU level into the subbasin level and then routed through the stream network.The Manning equation is used to estimate the flow rate and velocity through the channels.Flow routing is based on either the variable storage or the Muskingum routing method; and for this study, we chose the variable storage method (Neitsch et al., 2011).

Climate data
Hourly climate data for the period from 1998 to 2011 were measured and collected from several regional stations of the Korean Meteorological Agency (KMA) (Fig. 1).Precipitation and minimum/maximum temperature were obtained from the Haean KMA station (38.287 • N, 128.148 • E).Relative humidity, temperature, and wind speed were obtained from the Inje KMA station in the adjacent Yanggu County (38.207 • N, 128.017 • E).Solar radiation was collected from the Chuncheon KMA station (37.904 • N, 127.749 • E).Distributed climate data were also collected from 15 micrometeorological stations (Delta-T Devices, Ltd.) throughout the catchment (Fig. 1) between 2009 and 2011.Sub-hourly data was aggregated into hourly precipitation (±0.2 mm), minimum/maximum air temperature (±0.2 • C), wind speed (±0.1 m s −1 ), solar radiation (±5 W m −2 ), and relative humidity (±2 %).Each parameter was quality controlled by removing erroneous data and then gap filling from a similar station using a weighted algorithm based on elevation, station proximity, and aspect.The algorithm, as formulated for precipitation, is presented as dx dx +dy • P x −P y + P x w 2 + ...
The variable P e is the estimated precipitation (mm), z is the elevation (m), d is the distance to the observation point (m), φ is the observation point aspect (deg.), i is the time step, a is the total number of consecutive missing data, P o is the observed precipitation (mm), v is the total number of observational meteorological stations, j is the cumulative number of stations, the "e" and "o" subscripts are the estimated and observed location values, w is the weighting factor, and x and y subscripts are the first and second most proximal locations to the estimation location, respectively.Locally based relative humidity was modified by accounting for the temperature dependent local dew point.The SWAT model does not explicitly interpolate spatial meteorological conditions but instead, prescribes the nearest weather station parameters to the centroid of each subbasin (Neitsch et al., 2011).Due to the large variation in topographical complexity throughout the catchment, precipitation volume, soil moisture, and plant growth were impacted when SWAT assigned the meteorological data to each subbasin.We tested several interpolation methods to grid the measured meteorology results throughout the catchment (inverse distance weighted (IDW), spline, nearest neighbor, and kriging).The IDW method performed optimally and was used to grid the measured meteorological results throughout the catchment and the virtual weather corresponding to each subbasin centroid was prescribed.Principle data sources used for the Haean catchment ecohydrologic model are provided in Table 1.Choi et al. (2010) found highly variable temperature lapse rates, implying that stagnant East Asian monsoon high pressure systems can significantly vary climatic conditions on a local scale.A temperature lapse rate of −0.52 • C 100 m −1 was incorporated into the continuous spatial interpolation for temperature.

Discharge and evapotranspiration estimates
Event-based and baseflow surface water discharge measurements were collected at up to 14 locations throughout the catchment between 2003 and 2011 (Fig. 1) through multiple methods as described by Shope et al. (2013).Observed streamflow at interior locations within the catchment (S1, S4, S5, and S6) and the catchment outlet (S7) were utilized for daily and monthly model calibration to better parameterize spatial variability in hydrologic partitioning.These monitoring locations are distributed throughout the catchment along an elevation gradient with increasing drainage area and provide regional representation of model parameterization.In addition, the unique punchbowl shape enabled the calibration parameters to be correlated to other ungauged subcatchments with similar slope, elevation, and aspect.Spatiotemporal aquifer contributions were investigated by quantifying the relative baseflow from the hydrograph using several baseflow separation techniques including differential discharge measurements and recession analysis (Shope et al., 2013).For estimate consistency between each of the monitoring locations, we applied a recursive digital filter to separate the low-frequency baseflow signal from the high-frequency runoff in the formulation described by Eckhardt (2005).The calculated baseflow was subsequently compared to the SWAT modeled baseflow contribution.
The SWAT model also includes several methods to calculate potential evapotranspiration (PET) (Hargreaves and Samani, 1985;Monteith, 1965;Penman, 1948;Priestley and Taylor, 1972) depending on the observational meteorological data available.Because of the robust and highfrequency spatially variable micrometeorologic data available through the TERRECO project, we simulated daily PET using the Penman-Monteith method (Penman, 1948).As described in Ruidisch et al. (2013) and Shope et al. (2013), the weather conditions throughout the catchment are heterogeneous and therefore, the physically based Penman-Monteith estimates were preferred over the alternative methods.Soil evaporation and crop transpiration were estimated using the FAO Penman-Monteith equation as described in Allen et al. (1988).

DEM
The Soyang watershed 30 m resolution digital elevation model (DEM) obtained from the National Geographic Information Institute (NGII) was clipped to the extent of the Haean catchment boundaries (Fig. 1).The Haean catchment was divided into three slope classes representing steep forested high elevation (10 • to 90 • ), moderately sloped dry-land agriculture (2 • to 10 • ), and mildly sloping rice paddies in the central portion of the catchment (0 • to 2 • ) (Table 2).The observed river network was geo-referenced and explicitly incorporated into the DEM because modification of stream channels in highly managed catchments is prevalent and inclusion of stream delineation improves hydrologic segmentation and boundary delineation.In addition, extensive ground-based surveys of engineered channels, diversions, culverts, drainage features, sediment retention ponds, and roads throughout the Haean catchment were completed.To investigate the role that engineered structures have in channel routing, three channel classifications were constructed for (1) the river network; (2) the river network and engineered culverts; and (3) the river network, culvert system, and existing roads (Fig. 2).We implemented the engineered structures in SWAT by sequentially adding them to the prescribed river network and we superimpose the modified networks onto the DEM.The roads and culverts were then prescribed as impervious channels with no transmission loss on the river network.Therefore, we had three complete model constructs from the beginning to the end with different hydrographic segmentation and subbasin boundary delineation.

Soils
Regional soil information was obtained from the Rural Development Administration (RDA) (1 : 25000) and based on a single surficial soil layer.The Haean spatial soil data set  (TERRECO) coupled the RDA soil data, LULC, and extensive field-based soil profiles to develop a spatial distribution of multiple soil horizons to a depth of 3 m.Our results found that Haean soils are intensively managed and modified and highly dependent on land use (Tenhunen et al., 2011).Soil properties, including the hydrologic soil group, texture class, the percentage content of rock, sand, silt, and clay content, and the hydraulic conductivity, were derived from a 2009 catchment-wide field survey that was aggregated into 6 unique soil types (Table 2).The hydrologic group and texture for each of the soils is ( 1

Land use and land cover (LULC)
Intensive  a CN is the SCS curve number.b PHU is the cumulative heat units above 0.0 • C required for the LULC/crop to reach maturity.c Fertilizer type is classified as Chem (inorganic chemical) not explicitly described or Org (organic manure).d Fertilizer amount (kg ha −1 ).e Leaf out and cessation define the beginning and end of season for forest and orchard land use.

Management parameter estimation
Agricultural management practices within the Haean catchment were surveyed between 2009 through 2011 through a combination of on-site stakeholder interviews, empirical field observations (Tenhunen et al., 2011), published literature (i.e., Nguyen et al., 2012), and regulatory reports from the Research Institute of Gangwon (RIG), the Ministry of Environment, the National Institute of Agricultural Science, and Technology and the Korean Forest Research Institute.More than 300 interviews of stakeholders and farmers were completed under the TERRECO project to quantify fertilization and pesticide application quantities and timing, irrigation practices, planting and harvesting activities, and tillage methodologies.TERRECO managed plots were also used to obtain comprehensive temperature-based planting, fertilizer, tillage, mulching, development, and harvest information (J.Tenhunen, unpublished data).An example of the land use and crop management schedule, application rate, and appli-cation frequency is provided in Table 3. Fertilizer application parameters within the SWAT database were varied for each crop and subbasin for spatially distributed management.The simulated timing of management actions (i.e., fertilization, tillage, planting, irrigation, harvesting) was implemented in SWAT through daily heat unit summations because traditional planting and harvest methods are dependent on climatic observations closely correlated to heat units.

Biomass sampling, analysis, and plant growth
Biomass analysis was completed by collecting and sampling 5 to 10 entire plants, representative of each crop type (Table 2) from a 2009 catchment-wide sample set of TER-RECO harvest plots (J.Tenhunen, unpublished data).Each of the plants was field separated and subsequently weighed for fresh weight.The leaf area was individually measured using a portable leaf area meter (Opti-Sciences, Inc., AM 300).The samples were then separated and dried at 80  a Heat Units is the total base zero annual heat units for the plant cover/land use to reach maturity.b HUSC is the fraction of the total base zero annual heat units at which the management operation occurs.c BLAI is the maximum potential leaf area index.d DLAI is the fraction of the growing season when the leaf area begins to decline.e FRGRW1,2 represent the fraction of the plant growing season corresponding to the 1st and 2nd point on the optimal leaf area development curve.f LAIMX 1,2 represent the fraction of the maximum leaf area index corresponding to the 1st and 2nd point on the optimal leaf area development curve.g GSI is the maxixmum stomatal conductance at high solar radiation and low vapor pressure deficit (m s −1 ).h T_BASE is the minimum or base temperature for plant growth ( • C). i ALAI_MIN is the minimum leaf area index for the plant during the dormant period (m 2 m −2 ).j HVSTI is the fraction of aboveground biomass removed during a harvest operation and lost from the system.k CHTMX is the maximum canopy height (m).l BIO_E is the radiation use efficiency or biomass energy ratio ((kg ha −1 )/(MJ m −2 )).m BIO_LEAF is the fraction of tree biomass accumulated each year that is converted to residue during dormancy.n BMDIEOFF is the biomass die-off fraction.
To differentiate between crop types particular to South Korea (i.e., ginseng), several modified land use classes were created in the SWAT crop database.Nine representative field plots along an elevation transect were analyzed and crop parameters were varied to minimize the simulated and observed residuals for leaf area index (LAI), biomass, and crop yield.The crop parameters were altered based on observed measurements, plant physiology modeling results from the PIX-GRO model (i.e., Adiku et al., 2006), and published literature.The crop parameters that were varied are presented in Tables 3 and 4. Intensive cultivation was also present in agricultural areas not serviced by irrigation canals and therefore, groundwater abstraction was estimated from the PIXGRO model as the quantity required for optimal plant growth.Typical to many Asian catchments, Haean can be considered a highly managed catchment with increased uncertainty due to insufficient spatiotemporal water management data.

Rice paddies, potholes, and water abstraction
The quantity and timing of river and groundwater abstractions is uncontrolled and local estimates were inadequate for model inclusion.Depending on the HRU location, irrigation water was extracted from an adjacent river reach or from shallow groundwater.Groundwater-derived irrigation practices were limited to orchards and rice paddies and were accounted for in the simulations through water availability based auto-irrigation at the HRU level and defined by the soil water deficit.Haean rice paddies were simulated in SWAT as potholes, which are hydrologically similar to ponded areas.Rice paddies are typically characterized by multiple cascading-elevation plots separated by embankments.The rice paddies had low infiltration and typically saturated soil conditions and therefore, infiltration as a function of water content rather than flow routing was used for estimation of subsurface losses.The HRUs within each subbasin were developed using 0 % land use and 0 % slope threshold for reach subbasins resulting in maximum number of HRUs.Since a subbasin can have multiple HRUs but only have a single pothole, we limited the rice paddies in each subbasin to a single HRU.We accomplished this by varying the soil threshold until only a single rice paddy HRU was in each of the subbasins.

Meteorological drivers and the effects of interpolation
Meteorological time series data, particularly precipitation is a highly sensitive driver in hydrologic modeling applications (Strauch et al., 2012).Spatial monitoring distributions are typically limited and do not capture heterogeneous meteorological conditions that can be interpolated by wide-meshed monitoring networks (Notter et al., 2007).Large variations in elevation throughout the Haean catchment influence the precipitation volume, soil moisture, and plant growth.They can also influence the peak flow and the time of concentration to peak discharge of the simulated hydrograph (Khakbaz et al., 2012;Wilson et al., 1979).Our weather analysis revealed heterogeneous meteorological conditions throughout the Haean catchment that are dependent on elevation and aspect and largely focused in subregions (Choi et al., 2010;Shope et al., 2013).These meteorological variations have a direct influence on the relative humidity and therefore, the spatial variability of plant growth parameters between subbasins was significant (Fig. 3).
We examined the model sensitivity to alternative precipitation interpolation methods (IDW, Spline, nearest neighbor, and kriging), both through spatially explicit plant growth response and river discharge to assess the robustness of interpolation in our domain.We found that total river discharge between interpolation methods varied less than 0.1 % at the integrated catchment outlet (S7) and the discharge differences at multiple locations throughout the catchment (S1, S4, S5, and S6) were negligible.The IDW univariate interpolation technique for precipitation did result in slightly  improved plant growth response for selected crops and locations than other methods.Similar to results obtained by Notter et al. (2012), the IDW method was invoked to develop a continuous grid of meteorological drivers that were subsequently assigned to individual subbasins.

Sensitivity and model parameterization
The model sensitivity was addressed with respect to spatial distribution (number and location of meteorological stations, LULC distribution), observational record (LULC coverages, meteorological stations), resolution (soil coverage, subbasin discretization), and hydrologic stimulus (rainfallrunoff).The Haean catchment model configuration resulted in 142 topographically based subbasins and 2532 individual HRUs.Previous investigations have shown that the number of subbasins has little influence on runoff (Jha et al., 2004;Tripathi et al., 2006;Xu et al., 2012a, b).Alternatively, other studies have found that HRU discretization can have a substantial effect depending on the physical catchment conditions, data quality, and investigative scale (i.e., Setegn et al., 2008;Haverkamp et al., 2002).We assessed the effect of subbasin size and HRU definition on surface water dis-charge and found no appreciable difference between model results.However, our results show that elevation-based plant parameters and convective precipitation captured through increased subbasin discretization can be important.Subbasins with steep slopes and extensive vertical gradients must account for elevation-based climate conditions, which contribute to highly variable ET conditions.The sensitivity analysis of discharge related model parameters was achieved by sequentially varying an individual parameter while maintaining the remaining parameters for each monitoring location.
Between eight and eleven parameters from the original 15 discharge-related parameters were found to be sensitive to catchment-wide flow partitioning (Fig. 4).Subsequently, the range of each of the parameters was minimized during calibration procedures.The use of lumped, semi-distributed, and fully distributed model parameterization was also investigated through sensitivity analysis.We assigned the same parameter magnitudes by crop type for the lumped distributed parameters, by crop type and subbasin for semi-distributed, and by HRU in the fully distributed construction.We found that fully distributed parameters between subbasin, soil, and LULC were negligibly better than semi-distributed parameters based on aggregated LULC within individual subbasins.We also found that the use of a lumped parameter assignment did not perform as well as either the fully or semi-distributed parameterization.Therefore, for computational efficiency, a semi-distributed approach was taken throughout the catchment utilizing the most sensitive parameters at each monitoring location for parameterization in adjacent areas.
While we did not explicitly quantify the optimal parameterization, through a series of iterations we weighted the objective functions (R 2 , NSE, PBIAS, and baseflow percentage) in decreasing order as we compared individual locations throughout the catchment.In effect, we used a multi-criteria decision making process to determine the relative priority of each alternative when all of the criteria were considered simultaneously.Because our results indicated that the sensitivity analysis was significantly based on the monitoring location, we calibrated multiple locations along an elevation transect.In Fig. 4, the "t stat" provides a measure of parameter sensitivity where larger absolute values are more sensitive and the "p value" determines the significance of sensitivity with higher significance as values approach zero (Abbaspour, 2011).
Our results generally indicate surface runoff and routing parameters are more sensitive at higher elevations with increasing sensitivity to infiltration and groundwater parameters at lower elevations (Fig. 4).The REVAPMN groundwater parameter was a sensitive parameter at each location; however, the magnitude was relatively small.CH_K(2) was the least sensitive parameter, although included in the analysis for comparison.Table 6 provides a summary of the SWAT parameters.The infiltration parameters suggest significant baseflow response at higher elevations.At mid-elevations, surface runoff and routing parameters become more sensitive.At lower catchment elevations, infiltration, routing, and groundwater parameters dominate.Since the upper elevation locations are composed of shallow, highly permeable (S.Arnhold, unpublished results) soils over bedrock; we conceptualize high infiltration rates that contribute to increased baseflow and streamflow accumulation.At mid-to low-elevation locations, higher land management and soil amendments lead to runoff and less infiltration.These results identify the importance of and differences between model sensitivities as a function of the model equations, model sensitivity, and observational dynamics.Therefore, caution should be exercised in rainfall-runoff process simulations in relatively ungauged basins.

Metrics of model performance for calibration procedures
Model performance was assessed by several metrics at each location including the simulated and observed water balance, the coefficient of determination (R 2 ), Nash-Sutcliffe efficiency (NSE), percentage bias (PBIAS), and the baseflow contribution.The R 2 was used to evaluate time and space dependent cross-correlations and indicate if system behavior is accurately represented by the model (Bennett et al., 2012).The Nash-Sutcliffe efficiency (NSE) is a normalized correlation related statistic used to compare observational variance to the residual variance, particularly during peak events (Nash and Sutcliffe, 1970).The percentage bias (PBIAS) is a quantitative measure of simulated versus observed river discharge for the entire simulation period and defines the total volume differences between the simulated and observed fluxes.In addition, the baseflow statistic compares the simulated baseflow contribution to the calculated estimate at each location to alleviate hydrologic partitioning from alternative sources.This metric provides an independent check on a specific component of the water budget.Finally, measured plant growth dynamics were compared with simulated results.

Manual and automated model calibration
Due to the complexity of large-scale multi-objective analyses, watershed models are typically highly parameterized and manual calibration can be virtually impossible (Schuol and Abbaspour, 2006) although multi-site, multi-objective inverse calibration and uncertainty analysis can aid in understanding the system (Abbaspour et al., 2004;Duan et al., 2003).Model calibration was separated into two components, (1) manual catchment-scale calibration to estimate system processes and variability, and (2) automated calibration to quantify model uncertainty.
The SWAT model was simulated from 2006 through 2011 with the first 3 yr excluded for model initialization.The calibration and validation of river discharge was performed at a daily time step from 2009 through 2011, with 2010 as the calibration period and 2009 as the validation period.For locations S4 and S6, we did not have observational records for the 2009 validation period and instead used the concept of selfsimilarity for validation results.Since the transect followed an elevation gradient in a limited portion of the catchment, we conceptualized that similar hydrologic processes were occurring for similar elevation and drainage areas in other parts of the catchment.For example, location S4 was calibrated to the 2010 observational data, although there was limited data to validate for 2009.Because SD and SK had similar topography, elevation, drainage area, and land use patterning as S4 and S6, respectively, they were used to validate the S4 calibration parameters.Intensive manual calibration was performed at each of the subbasins routed to a monitoring station and used to minimize the acceptable parameter range at each site.The difficulty is that manual calibration sensitivity suffers from the linearity assumption by not accounting for correlations between individual parameters.
After manual calibration was optimized through the weighted, multi-criteria metrics previously discussed, automated model calibration, validation, and uncertainty analysis was completed using the Sequential Uncertainty Fitting Algorithm (SUFI-2) (Abbaspour et al., 2004(Abbaspour et al., , 2007)).The manual calibration results provided distributed, physically based parameter ranges that were incorporated into the SUFI-2 auto-calibration routine, starting with the catchment outlet and following a top to bottom approach.Model uncertainty in auto-calibration is quantified by the 95 % prediction uncertainty (95PPU) at the 2.5 and 97.5 % cumulative distribution, which is obtained through Latin hypercube sampling procedure (Abbaspour et al., 2004).Because the model varies multiple parameters at the same time, two indices are used to assess the stochastic calibration performance.The "p factor" describes the percentage of data bracketed by the 95 % prediction uncertainty and the "r factor" describes the Table 6.Calibration and validation statistics for each of the monitoring locations throughout the Haean Catchment.The data includes the subbasin demarcation of the monitoring locations, the total number of observations, the observed and simulated water balance, the NSE, R 2 , and PBIAS statistics, and the percent baseflow contribution.average width of the prediction band divided by the standard deviation of the measured data (Faramarzi et al., 2009).Since the uncertainty in field-based river discharge measurements was typically < 5 % (Shope et al., 2013), a conservative 10 % measurement error was included in the "p and r factor" calculations (Abbaspour et al., 2009;Andersson et al., 2009;Butts et al., 2004;Schuol et al., 2008).Yang et al. (2008) found that reasonable prediction uncertainty ranges were achieved with 1500 model simulation iterations, while, (Güngör and Göncü, 2012) showed that 300 iterations provided similar results to 1500 iterations.In Haean, at least 300 simulation iterations at each location were performed throughout the auto-calibration routine (Table 5).
As described, the calibration parameters were selected to optimize the PBIAS, R 2 , and NSE test statistics, the estimated groundwater baseflow, and the plant growth dynamics.The main SWAT parameters controlling baseflow processes in Haean include GW_REVAP, GWQMN, GW_DELAY, ALPHA_BF, and ESCO (Table 6).The primary parameters that affected surface runoff throughout the Haean catchment are CN2 and SOL_AWC.During model calibration procedures, the ESCO and GW_REVAP parameters were typically adjusted to minimize the PBIAS and improve the annual discharge and water balance trends.The GWQMN parameter was then adjusted to simulate the seasonal discharge trends assessed by maximizing the monthly R 2 and NSE statistics.Finally, the CN2, CH_N(2), and GWDELAY parameters were calibrated to account for daily trends by maximizing the NSE.When the Muskingum routing method was utilized, the channel parameters CH_N(2) and CH_K(2) were ranked 2 and 3 in the sensitivity analysis.However, the relative change in NSE between outlet results was negligible (∼ 0.01) compared to the default variable storage routing method, and the addition of more parameters was substantial.Therefore, variable storage routing within the SWAT model was chosen to limit the model parameterization.
The explanation for the deviations in runoff at the low elevation locations (S6 and S7) is not known or reflected in the SWAT input data.However, by examining a combination of optimized calibrated data, process-based comparisons, and field observations, the overall calibration metrics indicated increased flow routing directly from high elevation locations to lower elevation river locations.A possible explanation is the density of surface water collection and sedimentation ponds within the catchment, which may have impacted the observed runoff characteristics of the watershed (Cho et al., 2012).Using a multi-criteria optimization approach, we identified that engineered flow routing and infrastructure construction such as roads and culverts, contributed to increased discharge at lower elevations.These catchmentwide landscape engineering results are further discussed in Sect.4.5.

Spatiotemporal flow partitioning with respect to river discharge
The calibration and validation of the Haean catchment daily discharge yielded good results given the scarcity and the temporal longevity of the available data.The modeling results indicated that SWAT performance at the Haean catchment relied heavily on the quality and more importantly abundance of discharge data, similar to the results of Dessu and Melesse (2012).The NSE score for monitoring locations S1, S4, S5, S6, and S7 ranged between 0.64 and 0.95 with an average score of 0.76 for the 2010 calibration period and between 0.40 and 0.98 for the validation period (Fig. 5).Satisfactory NSE scores of > 0.5 (Moriasi et al., 2007) were achieved at all 14 gauge locations in the calibration period and at 12 of 14 in the validation period.The R 2 value was also reasonable for each of the monitoring locations, ranging from 0.70 to 0.96 with an average value of 0.81 for the calibration period and between 0.71 and 0.97 for the validation period (Fig. 5).
The fact that similar performance measures were reached in both validation and calibration periods indicate that there was minimal "overfitting" of the distributed parameters.
The baseflow contribution estimated at monitoring location S4 using a digital filter hydrograph separation technique was 26 %, although the calibrated estimate was 42 %.The hydrograph separation magnitude varied significantly, depending on the data quality, the length of the analysis, and the time step investigated.However, the digital filter methodology for estimation of the hydrograph separation is not process-based and may have significant uncertainty.The calibrated baseflow of 42 % at S4 is similar to the estimate at the upstream location S1 and nearly twice as high as all of the downstream locations, indicating that this mid-elevation area may be transition zone between baseflow and runoff dominated streamflow.This suggests that high elevation locations have increased baseflow contributions, relative to low elevation locations, regardless of the observational data period.
We found increased differences between the simulated and observed water balance as measured through PBIAS statistics at locations S6 and S7, which were 41 and 29 %, respectively.These PBIAS estimates are unsatisfactory according to Morasi et al. (2007), regardless of the very good R 2 and NSE metrics and acceptable baseflow estimates.The increase in water balance was hypothesized to be a function of rapid and large flow contributions from high elevation locations that were routed through culverts, drainages, and road systems to lower catchment locations.Essentially, the effect of the anthropogenic routing not only creates a large disparity in simulated discharge, but limits the subsurface infiltration at the plot-scale for higher elevation locations and surreptitiously develops a misleading flashy flow system with reduced landscape water storage.
The lower NSE score and R 2 values could be attributed to the low magnitude relative variability of discharge at higher elevation monitoring locations, which contributes to increased deviations of NSE scores during event conditions, particularly monsoonal extreme events.At location SK, there is scarce observation data and because the NSE statistics weight extreme events higher, limited but high deviations have a much larger impact than minor deviations.In addition, the difficulty in accurately simulating the river discharge at monitoring location SK was hypothesized to be a function of high elevation flow contributions that bypassed the monitoring gage as hyporheic flow (Shope et al., 2013).The hydrological response throughout East Asia and within the Haean catchment in particular, is typically flashy and erratic, further attributing to event-based deviations in the objective functions.At monitoring location S5, a higher temporal density of observations was obtained and the model performance metrics are generally better than for other locations.
Overall, the calibration and validation results were good and the percentage of baseflow contribution at each location was reasonable in terms of the hydrograph separation estimates.The auto calibration metrics of p-value and r value are both reasonable, while the R 2 and NSE statistics were consistently above satisfactory and predominately considered very good.The average p factor throughout the calibration period at all stations was 0.64 (0.54 to 0.69) and the r factor was 0.21 (0.10 to 0.38).The average p factor and r factor from the validation period was 0.74 (0.64 to 0.79) and 0.14 (0.10 to 0.21), respectively (Fig. 5).This indicates that the majority of the simulated results were within the 95 % confidence interval and that the standard deviation was adequately minimized.As shown in Fig. 5 and Table 7, the validation results  at these locations were good and consistent with the results estimated at the calibration locations.

51
Each of the objective functions, hydrologic partitioning quantified by PBIAS, and the baseflow percentages were calibrated simultaneously, which while optimizing the values of some parameters, were at the detriment of other parameters.For example, the NSE at S5 was initially 0.89; however, parameter adjustments were made to minimize the water balance, which resulted in a lower NSE value.The event on 1 September 2010 had a major influence on the magnitude of the NSE and R 2 objective function.This is primarily due to the paucity of observation points and therefore, the weight of individual points on the overall relationship, particularly during peak events.
The simulation results were very good in terms of adequately simulating baseflow contributions, the majority of moderate events, and most extreme events for each location.In addition, the other statistical objective functions were typically good to very good.The quality of input data, such as the estimated river discharge (Shope et al., 2013) or the short duration of observational data, significantly affected the model performance.For example, extensive observational data was collected at S5 but more limited at S4 and S6, resulting in decreased statistics at the latter location, even after calibration.The relatively large 95 PPU band "r factor" necessary to bracket the observed data indicates that the uncertainty in the conceptual model is also very important for the Haean catchment.

Agricultural management and production
The heat sum methodology used to estimate time variable management and planting actions, provides the flexibility to account for unseasonable variations in meteorological drivers between years (Fig. 3).Heat sums are calculated as the cumulative daily temperature greater than the base temperature of 0.0 • C initiated on the planting date and completed at the maximum growth.The HUSC is the percentage of the total heat units necessary for optimal growth of an individual crop and is prescribed for each management activity.The minimum heat sum over the period of record was 4246 • C during 2009, the maximum was 5783 • C during 2003, and the average annual heat sum is 5222 • C (Fig. 6).The 12 yr linear trend line of maximum cumulative annual heat sum values indicates a general decrease of nearly 74.8 • C per year.When the potentially extreme years of 2003, 2008, and 2009 were excluded, a decrease of 15.3 • C per year was estimated.While precipitation trends suggest more extreme events occurring over a shorter time, these results indicate a decreasing trend in annual heat output necessary for optimal plant growth.
To evaluate the SWAT simulation results on the ecohydrologic response, we also analyzed the simulation results in terms of agricultural growth dynamics at selected plot locations throughout the catchment.While calibrating spatiotemporal discharge as previously described, we also investigated the effect of crop dynamics through temporal leaf area index (LAI) as a proxy for crop growth and development (Fig. 7).Individual crop growth and development parameters were adjusted for a comparison between observed and simulated LAI (Table 4).Results indicate a generally reasonable approximation of simulated LAI where the R 2 for each of the crop types ranged from 0.51 to 0.76 (Fig. 7).More importantly, the results provide a consistent estimate of temporal trends in simulated biomass or agricultural production.

Influence of engineered landscape structure
Both the calibration and validation indicate successful spatial results with very good metrics, although a point of concern between observed and simulated results was at monitoring locations S6 and S7.The river discharge discrepancies between simulated and observed results were realized through PBIAS, which accounts for observed and simulated water balance differences.Field-based observations showed that catchment-wide surface runoff near the high elevation crops is routed to culverts immediately adjacent to the individual fields and road networks that discharge to low elevation river network reaches.As indicated in Fig. 2, many of these long, extensive features traverse from high elevation plots near the forest boundary down to the lower portions of the catchment.To test the impact of these anthropogenic engineered structures on catchment-wide hydrologic partitioning, we compared several different surficial flow routing configurations.The routing configurations utilized in the model simulations were (1) with rivers only, (2) with both rivers and culverts, and (3) a combination of rivers, culverts, and roads C. L. Shope et al.: Landscape complexity and ecosystem modeling with the SWAT model (Fig. 2).As previously described in Sect.4.3, the model performance in terms of PBIAS decreased toward the catchment outlet, particularly near S6 and S7.As the transect continues to the catchment outlet, the p factor decreases from 71 to 11 %, indicating that less data is bracketed by the 95 % confidence interval, while the r factor describing the standard deviation of the observed discharge increases from 0.20 to 0.36.
When the model was reconfigured to account for both the river drainage network and the culverts, a better calibration was obtained where the PBIAS at monitoring locations S6 and S7 decreased from 41 and 29 % to 8 and 9 %, respectively.The dramatic difference in PBIAS was not extended by including the roads into the river and culvert drainage network with a negligible increase in PBIAS observed at S6 and S7.Therefore, inclusion of the field-based drainage culverts was effective in moderating the difference in observed and model computed river discharge at lower elevation monitoring points and consistent with field-based observations of event-peak flow routing through the Haean watershed.However, it is surprising that the road network had minimal influence.During peak event conditions, substantial overland flow and sediment transport was observed throughout the Haean catchment.Since the poured concrete culverts are immediately adjacent to many of the plots, reduced landscape-scale infiltration required to maintain local soil moisture storage and rapidly transported excessive nutrients from fertilizer applications into the lower parts of the catchment is prevalent.This results in a rapid transport of elevated nutrient and sediment loads into the river.Therefore, while there is a significant influence on landscape-scale surface runoff, river discharge, and effectively hydrologic partitioning, a potentially greater issue is the impact expected from the rapid and large-scale alteration in water quality.

Conclusions
To provide a high accuracy estimate of spatiotemporal meteorological conditions, we used a unique high-frequency, quality control, and gap-filling algorithm to develop a detailed interpolation of weather patterns.The interpolated meteorological conditions were then discretized throughout the catchment and the conditions were prescribed at the centroid of each of the subbasins.This novel technique provided a better estimate of the dynamic variability due to convective storm events than the default SWAT application of prescribing the nearest weather station to the subbasin centroid.
We demonstrate that the use of a novel catchment-wide, multi-location, multi-objective function approach can drastically improve process-based estimates of catchment-wide hydrologic partitioning.By calibrating the model to many locations distributed throughout the catchment, landscape controls on hydrologic partitioning can be estimated as opposed to the integrated effect simulated at the catchment outlet.Be-cause the catchment is essentially a bowl-shaped topographic feature, the concept of symmetry enabled the results from a single elevation-based transect of monitoring locations to be utilized in a catchment-wide model calibration and validation.Our results showed that a combination of statistical, hydrologic, and plant growth objective functions as modeling metrics provide a more comprehensive understanding of system interactions.We included not only classical statistical metrics to calibrate our model, but we also calibrated the model to independent baseflow contribution estimates and plant growth dynamics.These novel calibration metric additions enabled us to improve the simulated hydrologic partitioning distributed throughout the catchment.
Our goal of simulating high-frequency monsoonal events in an area of complex physiographic topography provided substantial reliability in the use of the SWAT model in similar mountainous areas, particularly throughout East Asia.To enhance the calibration of the SWAT model, simulation of daily spatiotemporal stream discharge was improved through the incorporation of additional modeling metrics.Spatial variations of baseflow contributions and spatiotemporal plant growth dynamics through LAI helped to better constrain catchment-wide hydrologic partitioning.Our results show that fundamental shifts between surficial and baseflow driven hydrologic flow partitioning occur within the catchment.High elevation steep sloping regions were found to be generally baseflow dominated while lower elevation locations were predominately influenced by surface runoff.
The influences of engineered infrastructure systems (roads and culverts) were significant in hydrologic flow partitioning.Our results indicate that multiple calibration metrics and hydrologic characteristics (R 2 , NSE, PBIAS, baseflow percentage, and plant growth) were influential in quantifying scale-dependent watershed processes.By not including the culverts into the simulations, we demonstrate that the model simulations adequately represented observed spatiotemporal discharge.However, by including PBIAS as a calibration metric, we improved flow partitioning on the landscape scale by up to 33 %, particularly at the low elevation locations while minimal variations were observed at upper elevations.To optimize PBIAS, we explicitly included the culverts and the culverts and roads into the modeled drainage system to demonstrate that the spatially extensive irrigation culverts adjacent to most fields and the road network play an important role in flow routing.
However, there were limitations in the reliability of modeling in similar regions, particularly with respect to field estimates, data collection, and the conceptual model.In relatively ungauged locations, it can be difficult to adequately distribute a monitoring network with high-frequency temporal resolution.Data gaps due to equipment malfunction and instrument sensitivity to ice can be prevalent in locations with complex topography and meteorological variability.Another significant source of uncertainty is irrigation and consumptive use water withdrawal quantification.However, limited detailed data is typically available on the quantity, timing, or location of water withdraws and care should be taken to incorporate into model construction.
Overall, the results of this study show that unique modeling methodologies can be employed to decrease modeling uncertainty including accurate meteorological boundary conditions, spatially distributed monitoring locations, and additional physically based modeling metrics.Our results further elucidate the effect of catchment-scale engineered structures on discharge and the potential influence on nutrient loading and contaminant transport.Care must be taken during model construction to avoid overlooking valuable hydrologic information and complex relationships that may be deciphered through additional objective function metrics.This study shows the challenges of applying the SWAT model to complex terrain and meteorological extreme environments and the means to overcome these difficulties.

Fig. 1 .
Fig. 1.Haean study area within the Lake Soyang watershed is located in northeastern South Korea along the demilitarized zone (DMZ) border with North Korea.The regional KMA weather station and local meteorological stations are denoted with white circles and (WS).River discharge monitoring locations are denoted by (S) and the yellow squares.
Meteorologic variability and average daily value of each variable throughout the 1079 Haean catchment for 2010.Panel A) describes the daily precipitation and temperature 1080 variability, B) is the range in solar radiation and the average value between all of the 1081 locations, C) is the wind speed variability, and D) is the relative humidity range.1082

Fig. 3 .
Fig. 3. Meteorologic variability and average daily value of each variable throughout the Haean catchment for 2010.(A) describes the daily precipitation and temperature variability, (B) is the range in solar radiation and the average value between all of the locations, (C) is the wind speed variability, and (D) is the relative humidity range.

Fig. 5 .
Fig. 5. Calibrated and validated daily comparison of drainage area normalized observed and simulated river discharge along the elevation transect of monitoring locations S1, S4, S5, S6, and the catchment outlet S7.Included on each panel are the objective function and optimization statistics.

Figure 5 -Figure 6 -
Figure 5 -Calibrated and validated daily comparison of drainage area normalized observed and 1088 simulated river discharge along the elevation transect of monitoring locations S1, S4, S5, S6, 1089 and the catchment outlet S7.Included on each panel are the objective function and 1090 optimization statistics.1091

Fig. 7 .
Fig. 7. Comparison of simulated versus observed leaf area index (LAI) for five of the primary crops grown in Haean and the deciduous forest.

Table 1 .
Principle input data sets for the construction of the Haean catchment SWAT model.
k Farmer, county, administrative interviews and field-based plots

Table 2 .
Percentage of Haean catchment associated with the individual aggregated land use, soil, and slope classifications.The slope classification generally defines the difference between forest, dryland farming, and rice paddy systems throughout Haean.

Table 3 .
Agricultural crop management schedule including planting and harvest dates, fertilization dates, amounts, and type of fertilizer, tilling dates and method, SCS curve number for each crop, and the heat units required to reach maturity.

Table 4 .
Example SWAT model crop parameter database variations in the Haean model.

Table 5 .
SWAT parameter sensitivity and significance between discharge parameters throughout the Haean catchment (Fig.4).Calibrated SWAT parameters for the Haean catchment, including the individual ranking along the elevation-based transect, the minimum and maximum parameter values for all subbasins accounted for by each monitoring location, and the average calibrated parameter value.Because of the distributed nature of the Haean model, individual parameters varied depending on crop type, elevation, aspect and therefore, a specific parameter value is not available.

Table 7 .
Biomass production and crop yield statistics for South Korea and specifically, for the Haean catchment.: Ministry for Food, Agriculture, Forestry & Fisheries (MIFAFF), Korea Rural Economic Institute, Korean Statistical Information Service (KOSIS), Korea Agro-Fisheries Trade Corp. (aT), Yanggu statistical year-book 2003-2011 from the Yanggu County Office, FAOSTAT 2008, World Bank 2009. Sources