Comparison of Calibration Approaches of the Soil and Water Assessment Tool (SWAT) Model in a Tropical Watershed

: Hydrologic models are indispensable tools for water resource planning and management. Accurate model predictions are critical for better water resource development and management decisions. Single-site model calibration and calibrating a watershed model at the watershed outlet are commonly adopted strategies. In the present study, for the ﬁrst time, a multi-site calibration for the Soil and Water Assessment Tool (SWAT) in the Kelani River Basin with a catchment area of about 2340 km 2 was carried out. The SWAT model was calibrated at ﬁve streamﬂow gauging stations, Deraniyagala, Kithulgala, Holombuwa, Glencourse, and Hanwella, with drainage areas of 183, 383, 155, 1463, and 1782 km 2 , respectively, using three distinct calibration strategies. These strategies were, utilizing (1) data from downstream and (2) data from upstream, both categorized here as single-site calibration, and (3) data from downstream and upstream (multi-site calibration). Considering the performance of the model during the calibration period, which was examined using the statistical indices R 2 and NSE, the model performance at Holombuwa was upgraded from “good” to “very good” with the multi-site calibration technique. Simultaneously, the PBIAS at Hanwella and Kithulgala improved from “unsatisfactory” to “satisfactory” and “satisfactory” to “good” model performance, while the RSR improved from “good” to “very good” model performance at Deraniyagala, indicating the innovative multi-site calibration approach demonstrated a signiﬁcant improvement in the results. Hence, this study will provide valuable insights for hydrological modelers to determine the most appropriate calibration strategy for their large-scale watersheds, considering the spatial variation of the watershed characteristics, thereby reducing the uncertainty in hydrologic predictions.


Introduction
Hydrologic models are one of the most effective and widely accessible tools in water resource engineering and management. The introduction of computer science and its attributes to water resources engineering has paved the way for many computational models for hydrological studies. The Soil and Water Assessment Tool (SWAT) [1], HBV-light [2], and Hydrologic Engineering Center-Hydrologic Modeling System (HEC-HMS) [3] are to name a few. To obtain accurate predictions from hydrologic models, these should be appropriately calibrated, representing the specific watershed characteristics. It is acknowledged The Kelani River Basin stretches for 145 km and starts from Adam's Peak Mountains in the island's central hills. The focus basin in the present work extends approximately 2230 km 2 (which is slightly lower than the total catchment area of the KRB) discharging to the Indian Ocean in the Kelaniya area, a suburb in the Colombo District. Upstream of the river basin is green with forest, whereas the land use in the downstream part mainly comprises rubber plantations, homestead gardens, and built-up areas. The soil in the study region inherits a loamy nature characterized by low infiltration rates [39]. Figure 1 presents the topography of KRB with the selected meteorological stations (rain gauges and temperature gauges) and flow gauges. The river basin can be divided topographically into upper and lower basins. The upper basin is located upstream of the Hanwella River gauging station and spans approximately 1740 km 2 , while the lower basin is located downstream of Hanwella and covers about 500 km 2 . The lower basin is heavily urbanized, whereas the upper basin is mainly covered with dense vegetation such as tea, rubber, coconut, and forest [40]. The river discharges vary from 800 to 1500 m 3 /s during monsoon seasons, whereas it falls to 20 to 25 m 3 /s in the dry season, depending on the operation of three reservoirs in the catchment [41].
available in the watershed. Three scenarios of calibration were tested in this study. They include using (1) data from downstream, (2) data from upstream, and (3) data from downstream and upstream (multi-site calibration). Therefore, this study's primary objective is to comprehensively understand calibration methodologies and their effectiveness when calibrating small to large river basins.

Study Area
The Kelani River Basin stretches for 145 km and starts from Adam's Peak Mountains in the island's central hills. The focus basin in the present work extends approximately 2230 km 2 (which is slightly lower than the total catchment area of the KRB) discharging to the Indian Ocean in the Kelaniya area, a suburb in the Colombo District. Upstream of the river basin is green with forest, whereas the land use in the downstream part mainly comprises rubber plantations, homestead gardens, and built-up areas. The soil in the study region inherits a loamy nature characterized by low infiltration rates [39]. Figure 1 presents the topography of KRB with the selected meteorological stations (rain gauges and temperature gauges) and flow gauges. The river basin can be divided topographically into upper and lower basins. The upper basin is located upstream of the Hanwella River gauging station and spans approximately 1740 km 2 , while the lower basin is located downstream of Hanwella and covers about 500 km 2 . The lower basin is heavily urbanized, whereas the upper basin is mainly covered with dense vegetation such as tea, rubber, coconut, and forest [40]. The river discharges vary from 800 to 1500 m 3 /s during monsoon seasons, whereas it falls to 20 to 25 m 3 /s in the dry season, depending on the operation of three reservoirs in the catchment [41].

Soil and Water Assessment Tool (SWAT) Model
The Soil and Water Assessment Tool (SWAT) is a river basin model which was developed by the Agricultural Research Service (ARS) of the United States Department of Agriculture (USDA). The model can be used to assess the temporal and spatial behavior of flow situations in rivers. SWAT was initially used to analyze the impact of land uses and land management strategies on water, sediment, and chemical yields in complex watersheds [42]. It helped to work on extended and long-term analysis. SWAT is widely used

Soil and Water Assessment Tool (SWAT) Model
The Soil and Water Assessment Tool (SWAT) is a river basin model which was developed by the Agricultural Research Service (ARS) of the United States Department of Agriculture (USDA). The model can be used to assess the temporal and spatial behavior of flow situations in rivers. SWAT was initially used to analyze the impact of land uses and land management strategies on water, sediment, and chemical yields in complex watersheds [42]. It helped to work on extended and long-term analysis. SWAT is widely used globally and is considered an adaptable model that can incorporate a variety of environmental processes to enhance more efficient watershed management and the formation of more informed policy decisions [43]. SWAT divides a watershed into sub-watersheds, each further subdivided into hydrologic response units (HRUs) comprised of land use, land management, and soil characteristics [42]. Equation (A1) (in Appendix A) presents the mathematical formulation of SWAT which is based on the water balance. The Kelani River Basin's watershed boundaries, stream network, topography, and subbasins with sub-basin characteristics such as slope and slope length were defined using a digital elevation model (DEM) with a 30 m resolution. The resolution of the DEM is a crucial consideration when choosing a DEM for a SWAT model, as it has a considerable impact on the total length of the streamline, the slope of the main channel, the watershed area, and the delineation of the area slope. Buakhao and Kangrang [44] completed a comprehensive analysis of this topic in the application of SWAT modeling and more information can be found in the related reference.

Land-Use Land-Cover (LULC) Properties
The field observation data and high-resolution satellite images from Google Earth were used to perform the land-use land-cover (LULC) classification. Landsat images with a resolution of 30 m × 30 m were used for classification. The image (Satellite Name: Landsat 5, Sensor ID: TM, Acquisition Date: 7 February 1997 and Path/Row: 141/055) retrieved from the United States Geological Survey (USGS) Earth Explorer (https://earthexplorer.usgs. gov/ accessed on 5 September 2022) had a cloud cover of 8.0%. In order to categorize the land into five distinct categories, the land-use classification system developed by the United States Geological Survey (USGS) was used. More information on the USGS classification system can be found in Anderson et al. [45].
Semi-automated classification plugin in QGIS 3.10 was used for the classification. QGIS is considered one of the most common and globally accepted open-source GIS applications. In addition, the training areas and pixel-based image classification were used to perform supervised classification. More information on this open-source toolbox can be found in Congedo [46]. Then, standard training samples were used to classify each land-use class following Lillesand et al. [47]. In addition, the QGIS semi-automatic classification plugin was utilized to evaluate land-cover classification accuracy. Stratified random points (a new function in SCP 6.4.0) were used in the SCP function to generate Region of Interests (ROIs), which were then photo-interpreted and used as a reference for the accuracy assessment. The error matrix and estimated accuracy values were then developed by the SCP Accuracy tool.
There are a variety of approaches for determining the appropriate sample size and distribution in these types of studies. As a general rule, random sampling is the most effective method for ensuring low estimates of standard errors in accuracy. The proportions of land cover classes and the standard errors that were anticipated for overall land cover categorization and individual classes significantly impact the sample design. Class-specific estimations can be improved by stratifying the sample. Olofsson et al. [48] provided further information on sample size and stratification. The sample size (N) can be computed as per Equation (1).
where W i is the class i area to total area (proportion); S i is the standard deviation of i stratum; S 0 is the standard deviation expected for overall accuracy.

Soil Properties
The soil properties of KRB were derived from the Food and Agriculture Organization (FAO) and The United Nations Educational, Scientific and Cultural Organization (UNESCO) soil maps. Four major soil types, loam, clay, and sandy soils (refer to Table 1) can be identified as per the SWAT global soil database. The soil map for the KRB is given in Figure 2. The spatial variation of soil types can be clearly seen in this figure. The upstream of the KRB mostly has Ao73-2bc-3645 soil whereas the downstream has both clay and sandy loam soils.

Soil Properties
The soil properties of KRB were derived from the Food and Agriculture Organization (FAO) and The United Nations Educational, Scientific and Cultural Organization (UNESCO) soil maps. Four major soil types, loam, clay, and sandy soils (refer to Table  1) can be identified as per the SWAT global soil database. The soil map for the KRB is given in Figure 2. The spatial variation of soil types can be clearly seen in this figure. The upstream of the KRB mostly has Ao73-2bc-3645 soil whereas the downstream has both clay and sandy loam soils.

Meteorological and Hydrological Data
Daily resolution meteorological data of rainfall, minimum and maximum temperature, relative humidity, solar radiation, and wind speed data have to be fed to the SWAT model. These data can either be ground measured or generated by the SWAT WXGEN weather generator. This weather generator helps to fill the gaps in ground-measured data [49]. Sri Lanka has widespread rainfall measuring stations; however, it only has a few temperature and other meteorological data measuring stations (less than 25 for the whole country). Therefore, observed rainfall and temperature data were used in modeling; however, due to lack of measured data, the relative humidity, wind speed, and solar radiation data from 1994 to 2015 were generated. The daily rainfall data of Laxapana, Hanwella, Wewalthalawa, Norwood, Kitulgala, Holombuwa, Deraniyagala, Glencourse, and Angoda, and temperature data of Colombo, Katunayake, and Rathnapura were purchased from the Department of Meteorology, Sri Lanka, from 1994 to 2015. These gauging stations are shown in Figure 1. In addition, daily discharge data for five flow stations, including

Meteorological and Hydrological Data
Daily resolution meteorological data of rainfall, minimum and maximum temperature, relative humidity, solar radiation, and wind speed data have to be fed to the SWAT model. These data can either be ground measured or generated by the SWAT WXGEN weather generator. This weather generator helps to fill the gaps in ground-measured data [49]. Sri Lanka has widespread rainfall measuring stations; however, it only has a few temperature and other meteorological data measuring stations (less than 25 for the whole country). Therefore, observed rainfall and temperature data were used in modeling; however, due to lack of measured data, the relative humidity, wind speed, and solar radiation data from 1994 to 2015 were generated. The daily rainfall data of Laxapana, Hanwella, Wewalthalawa, Norwood, Kitulgala, Holombuwa, Deraniyagala, Glencourse, and Angoda, and temperature data of Colombo, Katunayake, and Rathnapura were purchased from the Department of Meteorology, Sri Lanka, from 1994 to 2015. These gauging stations are shown in Figure 1. In addition, daily discharge data for five flow stations, including Hanwella, Glencourse, Kitulgala, Holombuwa, and Deraniyagala, were purchased from the Department of Irrigation, Sri Lanka. These data were used to conduct the sensitivity analysis first and then to calibrated and validate the SWAT model. modules for in-depth analysis. Initially, the input raster files (DEM, soil map, and LULC map) were all projected into the same geographic coordinate system (i.e., UTM zone 44N for Sri Lanka). Generally, the watershed delineation procedure consists of establishing a DEM, defining streams, inlets, and outlets, and calculating sub-basin parameters. The in-stream definitions process was not defined correctly by the ArcSWAT model and was not delineated through the Glencourse flow station. To address the issue with stream definitions, this study used the SAGA tool combined with the open street map (OSM) plugin in QGIS 3.10. The stream was burned to ArcSWAT after the river network was derived by identifying the Strahler order. A catchment boundary with a total area of 2239.5 km 2 was recognized by the model.
This KRB was subdivided into seven sub-basins with manually added outlets. After that, the land area of each sub-basin was split into hydrological response units (HRUs). Fifty-three HRUs were identified. The slope map, soil layers, and LULC were imported into the project using ArcSWAT's HRU analysis tool. Thus, the land use and soil layers have completely enclosed the delineated watershed. In addition to land use and soils, SWAT HRU analyses have included classifications of HRUs by slope classes; therefore, the multiple slope option was selected. After reclassifying the LULC, slope, and soil maps to correspond with SWAT database values, these physical attributes were overlaid to define the HRUs. In this study, the threshold values of 15%, 5%, and 5% were specified for land use, slope, and soil, respectively, for analysis.

Parameter Selection
The SWAT flow simulation was performed with monthly frequency after processing the required inputs from 1994 to 2015. The initial three years (1994 to 1996) were allocated to warm up the model. The SWAT model was performed with the default SWAT parameters during the simulation. Relevant SWAT parameters, including the baseflow alpha-factor (ALPHA BF.gw), the minimum water level required for the formation of baseflow (GWQMN.gw), the time needed for water to reach the shallow aquifer after leaving the root zone (GW DELAY.gw), the groundwater re-evaporation coefficient (GW_REVAP), the curve number (CN2.mgt), the threshold depth at which water can percolate from the shallow aquifer into the deeper aquifer (REVAPMN), the compensation factor for the soil evaporation (ESCO.hru), and the maximum canopy water storage (CANMX.hru) were selected by the literature and visual observations of both observed and simulated hydrographs. Overestimations, underestimations, and deviations in the discharge graphs were highlighted.
Before the manual calibration, a SWAT-CUP analysis was carried out to identify the most sensitive parameters. SWAT-CUP is an independent program designed for SWAT calibration [50]. Five distinct calibration processes (Particle Swarm Optimization (PSO), Markov Chain Monte Carlo (MCMC), Sequential Uncertainty Fitting ver.2 (SUFI-2), Generalized Likelihood Uncertainty Estimation (GLUE), and Parameter Solution (Parasol)) are included in the program, as well as functions for validation and sensitivity analysis and a Bing map visualization of the study area. As part of the SWAT-CUP analysis, this study adopted the Latin Hypercube Sampling technique from the SUFI-2 to calibrate and estimate uncertainty [50,51]. The SUFI-2 algorithm is the most common and computationally efficient method, having the best P-factor (prediction uncertainty ranges) and R-factor (relative coverage of measurements) compared to other methods [38,[52][53][54].
Even though there were seven outlets in the basin, only five were considered for the analysis. The fourth outlet, Nagalagam Street, is the closest outlet to the sea and is therefore subject to tides [54]. Consequently, it was not possible to observe or simulate the discharge, and outlet number 4 cannot be used (fourth sub-basin). In addition, the discharge data from Norwood (outlet 7) were inadequate for the analysis. The remaining outlets (refer to Figure 3) were calibrated and validated. The research included two calibration methods, single-site calibration and multi-site calibration.
Even though there were seven outlets in the basin, only five were considered for the analysis. The fourth outlet, Nagalagam Street, is the closest outlet to the sea and is therefore subject to tides [54]. Consequently, it was not possible to observe or simulate the discharge, and outlet number 4 cannot be used (fourth sub-basin). In addition, the discharge data from Norwood (outlet 7) were inadequate for the analysis. The remaining outlets (refer to Figure 3) were calibrated and validated. The research included two calibration methods, single-site calibration and multi-site calibration.
In single-site calibration, the authors used the upstream station (Deraniyagala) and the downstream station (Hanwella) and calibrated each case independently. The resulting sub-basin parameter values were then assigned to the remaining four sub-basins. Multisite calibration was accomplished by calibrating from the upstream to the downstream river while holding identical values for each sub-basin.

Performance Evaluation of the Model
The performance of the model was evaluated using statistical indices such as the Nash-Sutcliffe model efficiency (NSE), the RMSE-observations standard deviation ratio (RSR), the percentage of bias (PBIAS), and the coefficient of determination ( ) (find the equations and definitions in Appendix B). The ranges between 0 and 1 and measures the model's ability to explain the variance of observed data. The greater the number, the lower the error variance of the model, and vice-versa for lower numbers. Generally, > In single-site calibration, the authors used the upstream station (Deraniyagala) and the downstream station (Hanwella) and calibrated each case independently. The resulting sub-basin parameter values were then assigned to the remaining four sub-basins. Multi-site calibration was accomplished by calibrating from the upstream to the downstream river while holding identical values for each sub-basin.

Performance Evaluation of the Model
The performance of the model was evaluated using statistical indices such as the Nash-Sutcliffe model efficiency (NSE), the RMSE-observations standard deviation ratio (RSR), the percentage of bias (PBIAS), and the coefficient of determination (R 2 ) (find the equations and definitions in Appendix B). The R 2 ranges between 0 and 1 and measures the model's ability to explain the variance of observed data. The greater the number, the lower the error variance of the model, and vice-versa for lower numbers. Generally, R 2 > 0.5 is considered the acceptable range for a model [55,56]. The mathematical formulation for R 2 is given in Equation (A2).
The Nash-Sutcliffe efficiency (NSE) is a normalized statistic that evaluates the amount of residual variance relative to the variance of observed data and reflects the fit between simulated and observed data plots [57]. While Servat and Dezetter [58] state that the NSE is the best objective function to reflect the overall fit of a hydrograph, Legates and McCabe Jr [59] recommend using the NSE. When it is set to 1, the output is optimal. In general, models with an NSE of at least 0.5 are considered to be in acceptable form. The mathematical formulation of the NSE is given in Equation (A3).
Based on the recommendations of Singh et al. [60], the RMSE-observations standard deviation ratio (RSR) is computed as the ratio between the RMSE and the standard deviation of measured data. Root-mean-squared error (RMSE) is one of the most often employed error-index statistics [60,61]. Error-index statistics and a scaling/normalization factor make it possible to apply the final statistic and reported values to various constituents. With an RSR value between 0 and ∞, which shows 0 RMSE or residual variation and thus perfect model simulation, the RSR spans from 0 to a significant positive value. This produces a low RSR and decreased RMSE, which in turn improves model simulation outcomes. Equation (A4) presents the mathematical formulations of the RMSE-observation standard deviation ratio.
The percent bias (PBIAS) measures the average tendency of simulated data to differ from their observed counterparts in terms of percentage (refer to Equation (A5)). For example, a model with a negative PBIAS indicates overestimation, whereas a positive PBIAS indicates underestimating. PBIAS is frequently used to measure water balance errors and has the capability to reveal poor model performance [62]. During dry seasons, PBIAS values for streamflow tend to vary more among different autocalibration methods than in wet seasons [62]. Table 2 presents the summary of performance evaluation values for monthly discharge simulations, which were recommended by many researchers.

Land-Use Land-Cover Classification
The classification accuracy was examined using the Kappa coefficient and was above 81%, while the overall classification accuracy was over 84%, suggesting that the classification was acceptable for the analysis. Forest land area (1875.8 km 2 ) was the highest in 1997, while the water bodies (16.9 km 2 ) acquired a small percentage of the area in the Kelani River Basin. However, agricultural lands shared a significant percentage while more settlements can be seen downstream of the KRB. Figure 4 showcases the land use land cover for KRB in 1997. Other assessed land cover is given in Table 3 with their areas.     Initially, the SWAT model was calibrated using downstream data (Hanwella sub-basin), and the obtained values for the selected parameters were then assigned to other basins as is customary with the single-site calibration technique. Figure 5 showcases the calibration and validation of observed and simulated streamflows. The streamflow gauges Hanwella, Deraniyagala, Glencorse, Holombuwa, and Kithulgala illustrate a good agreement between observed and simulated streamflow data, highlighting a good performance of the model, while slight overestimations are seen at Glencorse and Kithulgala stations and slight underestimation can be observed during some of the months of 1997, 1998, 2010, and 2015 at Kithulgala station (refer to Figure 5c,e).
Throughout the calibration period, the model performances for Hanwella, Deraniyagala, and Holombuwa stations were deemed "good" with the indices RSR and NSE, values of (0.6, 0.65), (0.52, 0.66), and (0.56, 0.72), respectively (refer to Table 4). In addition, the PBIAS stayed less than ±10 at Deraniyagala and Holombuwa stations, and the model performance was assessed as "very good". However, the PBIAS ≥ 25 at Hanwella station was classified as "unsatisfactory". Although Kithulgala station shows a "satisfactory" condition with a PBIAS value of (−20.35), the model performance at Glencorse and Kithulgala is rated as "unsatisfactory" for all the statistical indices except R 2 . Overestimation of the peak flows could be the reason for low NSE values at Glencorse and Kithulgala stations because of the sensitivity of NSE to higher values resulting from squared differences in its derivation. The R 2 values ranged from 0.61 to 0.82, suggesting "very good" (0.75-0.89) to "good" (0.50-0.74) performance of the model [66]. During the validation period, model performance was "satisfactory", with NSE values ranging from 0.51 to 0.62 for all stations except Glencorse (−0.17), which can be rated as "unsatisfactory". The PBIAS and RSR can be classified as "very good" and "good" for Deraniyagala  Table 5 showcases the model calibrated parameters for Hanwella sub-basin and these values can be effectively used in any future study. Throughout the calibration period, the model performances for Hanwella, Deraniyagala, and Holombuwa stations were deemed "good" with the indices RSR and NSE, values of (0.6, 0.65), (0.52, 0.66), and (0.56, 0.72), respectively (refer to Table 4). In addition, the PBIAS stayed less than ±10 at Deraniyagala and Holombuwa stations, and the model performance was assessed as "very good". However, the PBIAS ≥ 25 at Hanwella station was classified as "unsatisfactory". Although Kithulgala station shows a "satisfactory" condition with a PBIAS value of (−20.35), the model performance at Glencorse and Kithulgala is rated as "unsatisfactory" for all the statistical indices except . Overestimation of the peak flows could be the reason for low NSE values at Glencorse and Kithulgala stations because of the sensitivity of NSE to higher values resulting from squared differences in its derivation. The values ranged from 0.61 to 0.82, suggesting "very good" (0.75-0.89) to "good" (0.50-0.74) performance of the model [66]. During the validation period, model performance was "satisfactory", with NSE values ranging from 0.51 to 0.62 for all stations except Glencorse (−0.17), which can be rated as "unsatisfactory". The PBIAS and RSR can be classified as "very good" and "good" for Deraniyagala varies from 0.61 and 0.85, signifying "very good" (0.75-0.89) to "good" (0.50-0.74) model performance. Table 5 showcases the model calibrated parameters for Hanwella sub-basin and these values can be effectively used in any future study.    Figure 6 presents the calibration and validation results for the Deraniyagala sub-basin. Similar to the initial single-site calibration approach, the simulated flow at Hanwella, Deraniyagala, Glencorse, Holombuwa, and Kithulgala gauging stations showed a good fit with observed data, with a slight overestimation of the flow at Glencorse and Kithulgala for the majority of the years of the studied period (refer to Figure 6c,e).

Deraniyagala Sub-Basin
In the calibration period, the model performances at Deraniyagala and Holombuwa stations were graded as "good" based on the NSE and RSR indices (0.66, 0.56) and (0.69, 0.57), respectively. In addition, Hanwella (0.65, 0.64) was rated "good" for the NSE and "satisfactory" for the RSR, while Kithulgala (−0.07, 1.06) was deemed "unsatisfactory" for both indices (refer to Table 6). The PBIAS at Deraniyagala and Holombuwa stations did not surpass ±10; hence they were assessed as "very good." Hanwella and Kithulgala were rated as "unsatisfactory" and "satisfactory" for the PBIAS, whereas Glencorse was rated as "unsatisfactory" for all indices except R 2 , as is customary. R 2 values ranged from 0.56 to 0.81, suggesting "very good (0.75-0.89)" to "good (0.50-0.74)" performance of the model.
The model performance for the validation period was "unsatisfactory" for the NSE, PBIAS, and RSR indices at Hanwella (0.41, −51.41, 0.87) and Glencorse (−0.15, −41.99, 1.06) stations. However, the values of PBIAS for Deraniyagala (6.06), Holombuwa (0.06), and Kithulgala (−7.96) stations were less than ±10, which indicates a "very good" model performance. The NSE and RSR indices demonstrated "good" performance at the Deraniyagala (0.68, 0.55) and "satisfactory" performance at the Kithulgala (0.55, 0.65) stations. Holombuwa's model performance with RSR was deemed "good", whereas the NSE suggests a "satisfactory" performance. R 2 indicates that the model performance at all stations falls within the range of "very good" (0.75-0.89) to "good" (0.50-0.74). Detailed analyses of these results are given in Table 6. Similar to the Hanwella sub-basin, Table 7 presents the calibrated parameters for SWAT model for Deraniyagala sub-basin. These calibrated values can be used in any future study.  Figure 6 presents the calibration and validation results for the Deraniyagala sub-basin. Similar to the initial single-site calibration approach, the simulated flow at Hanwella, Deraniyagala, Glencorse, Holombuwa, and Kithulgala gauging stations showed a good fit with observed data, with a slight overestimation of the flow at Glencorse and Kithulgala for the majority of the years of the studied period (refer to Figure 6c,e).  In the calibration period, the model performances at Deraniyagala and Holombuwa stations were graded as "good" based on the NSE and RSR indices (0.66, 0.56) and (0.69, 0.57), respectively. In addition, Hanwella (0.65, 0.64) was rated "good" for the NSE and "satisfactory" for the RSR, while Kithulgala (−0.07, 1.06) was deemed "unsatisfactory" for both indices (refer to Table 6). The PBIAS at Deraniyagala and Holombuwa stations did not surpass ±10; hence they were assessed as "very good." Hanwella and Kithulgala were rated as "unsatisfactory" and "satisfactory" for the PBIAS, whereas Glencorse was rated as "unsatisfactory" for all indices except , as is customary. values ranged from 0.56 to 0.81, suggesting "very good (0.75-0.89)" to "good (0.50-0.74)" performance of the model.

Deraniyagala Sub-Basin
The model performance for the validation period was "unsatisfactory" for the NSE, PBIAS, and RSR indices at Hanwella    The model performance for the Hanwella was assessed as "very good" based on the R 2 for both calibration and validation periods (refer to Table 8). The NSE and RSR indices exhibited "good" performance of the model within the calibration period, whereas the PBIAS was "satisfactory". However, throughout the validation period, the PBIAS and RSR evaluated the model performance in Hanwella as "unsatisfactory", although the NSE rated it as "satisfactory". In addition, the Deraniyagala sub-basin outperformed all other sub-basins, with all model performance indices falling within the range of "very good" and "good" for both calibration and validation periods. Holomabuwa follows a similar trend, except that the NSE graded the model's performance as "satisfactory" during the validation period. Furthermore, the model performance of Glencorse was graded as "unsatisfactory" by all indices except R 2 , which indicated a "very good" and "good" model performance for both calibration and validation, respectively. In the calibration period, model performance for the Kithulgala sub-basin was assessed as "good" based on the the R 2 and PBIAS indices. However, during the validation period, model performance was rated as "good" for both indices except the PBIAS, which increased to "very good" performance. NSE and RSR ranked the model performance as "unsatisfactory" and "good" in calibration and validation periods, respectively. Table 9 presents the fitted values for SWAT calibrated values. As it was presented under Tables 5 and 7, these calibrated values can be effectively used by any future study for the same catchment. The model performance for the Hanwella was assessed as "very good" based on the for both calibration and validation periods (refer to Table 8). The NSE and RSR indices exhibited "good" performance of the model within the calibration period, whereas the PBIAS was "satisfactory". However, throughout the validation period, the PBIAS and RSR evaluated the model performance in Hanwella as "unsatisfactory", although the NSE rated it as "satisfactory". In addition, the Deraniyagala sub-basin outperformed all other sub-basins, with all model performance indices falling within the range of "very good" and "good" for both calibration and validation periods. Holomabuwa follows a similar trend, except that the NSE graded the model's performance as "satisfactory" during the validation period. Furthermore, the model performance of Glencorse was graded as "unsatisfactory" by all indices except , which indicated a "very good" and "good" model performance for both calibration and validation, respectively. In the calibration period, model performance for the Kithulgala sub-basin was assessed as "good" based on the the and PBIAS indices. However, during the validation period, model performance was rated as "good" for both indices except the PBIAS, which increased to "very good" performance. NSE and RSR ranked the model performance as "unsatisfactory" and "good" in calibration and validation periods, respectively.  Table 9 presents the fitted values for SWAT calibrated values. As it was presented under Tables 5 and 7, these calibrated values can be effectively used by any future study for the same catchment.

Comparative Analysis
The average annual basin values from a single-site and multi-site calibration were obtained for comparative analysis and shown in Table 10. Both single-site calibrations showcase similar range values for average annual values. However, some significant changes can also be seen (surface runoff, groundwater, etc.). This considerable difference may result from the several calibration procedures with varying parameter values. For instance, the GW DELAY.gw values (146.45 and 340) were altered during two calibration operations, as illustrated in Tables 5 and 7, which can immediately impact the simulated groundwater flow into the watershed's stream. In addition, the multi-site calibration values are approximately similar to calibration 1 (Hanwella gauging station) from singlesite calibration. They do not show much difference. Hanwella gauging station is far downstream compared to the Deraniyagala gauging station, and this could be a reason for having a slightly increased deviation from single-site calibration to multi-site calibration.
The hydrological features of a basin are represented by its annual mean values for both calibration methods. The ratio between precipitation (P) and potential evapotranspiration (E p ) reveals the moisture state, which reflects the climate condition of the region. By analyzing this ratio, the climate condition of the relevant study region may be identified, and it may assist the authorities in determining whether the results of this study would accurately reflect the actual setting of the basin. According to the National Oceanic and Atmospheric Administration (NOAA), the relevant climate conditions and the ratio between P and E p are shown in Table 11.  This study indicates a P/E p value of 2.4 corresponds to expected humid conditions, given that the Kelani River Basin is a tropical watershed. This information is important for validating the remaining hydrological features of the basin.
A few studies on SWAT models in Sri Lanka have been conducted [68][69][70]. Among them, Siriwadana and Rajapakse [37] were the first to document SWAT model results in the Kelani River Basin. The results of our study will guide water resource management authorities in decision making with optimal data availability. Thus, the research work presented herein creates the novelty of the work in the context not only in the KRB but also in a country blessed with more than 100 rivers that are spread all over the country: Sri Lanka.
The results of this analysis incorporate several uncertainties. According to Gunathilaka and Samarasinghe [54,71], the backwater effect from tidal affects up to the Hanwella gauge (which is in the upstream). This can alter the discharge of the river. Moreover, the SWAT model has its own uncertainties. Since the model simulates watershed processes through several assumptions, significant uncertainties are carried out through the model simulation result.

Summary and Conclusions
This study evaluates the SWAT model calibration using three distinct calibration procedures, including (1) single-site calibration with downstream data, (2) single-site calibration with upstream data, and (3) multi-site calibration using downstream and upstream data approach for the first time in the context of Sri Lanka. Calibration was performed using SWAT parameters, including the baseflow alpha-factor (ALPHA BF.gw), the minimum water level required for the formation of baseflow (GWQMN.gw), the time needed for water to reach the shallow aquifer after leaving the root zone (GW DELAY.gw), the groundwater "re-evaporation" coefficient (GW_REVAP), the curve number (CN2.mgt), the threshold water level in the shallow aquifer for "re-evaporation" to occur (REVAPMN), the maximum canopy water storage (CANMX.hru), and the soil evaporation compensation factor (ESCO.hru). These parameters were selected by the literature and visual observations of both simulated and observed hydrographs.
The performance of the model was evaluated using statistical indices, including coefficient of determination (R 2 ), the Nash-Sutcliffe model efficiency (NSE), the RMSE observed standard deviation ratio (RSR), and the percentage of bias (PBIAS). These are the most prevalent statistical indices in the literature for evaluating SWAT model performance.
The results indicate a satisfactory fit between these calibration approaches (single and multi-site calibration) and the observed data, with a notable improvement in the multisite calibration strategy compared to the other two methods. The study concluded that, although the multi-site calibration strategy requires a considerable amount of time, the accuracy of the results can be improved with this method, and it will allow the hydrological modelers to consider the spatial variation of the watershed characteristics and reduce the uncertainty in hydrologic predictions.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A. Water Balance Equation
where SW t i is the ultimate water content (mm) at the time t, SW 0 is the initial water content of the soil (mm), t is the duration of the simulation (days), R day i is the amount of precipitation recorded on the ith day (mm), Q sur f i is the amount of surface runoff recorded on the ith day (mm), E a i is the amount of evapotranspiration recorded on the ith day (mm), W seep i is the infiltration into the vadose zone on the ith day from the soil profile (mm), and Q gw i is the amount of baseflow recorded on the ith day (mm).

Appendix B. Statistical Indices Utilized for Model Performance Evaluation
Here, Q represents a variable (e.g., discharge) and s, m subscripts represent the simulated value and the measured value. Therefore, Q m is the mean measured data and Q s is the mean simulated data. In addition, i describes the ith measured or simulated data while n represents the total number of measured or simulated data.