Significance of Multi-Variable Model Calibration in Hydrological Simulations within Data-Scarce River Basins: A Case Study in the Dry-Zone of Sri Lanka

: Traditional hydrological model calibration using limitedly available streamflow data often becomes inadequate, particularly in dry climates, as the flow regimes may abruptly vary from arid conditions to devastating floods. Newly available remote-sensing-based datasets can be supplemented to overcome such inadequacies in hydrological simulations. To address this shortcoming, we use multi-variable-based calibration by setting up and calibrating a lumped-hydrological model using observed streamflow and remote-sensing-based soil moisture data from Soil Moisture Active Passive Level 4. The proposed method was piloted at the Maduru Oya River Basin, Sri Lanka, as a proof of concept. The relative contributions from streamflow and soil moisture were assessed and optimised via the Kling–Gupta Efficiency ( KGE ). The Generalized Reduced Gradient non-linear solver function was used to optimise the Tank Model parameters. The findings revealed satisfactory performance in streamflow simulations under single-variable model validation ( KGE of 0.85). Model performances were enhanced by incorporating soil moisture data ( KGE of 0.89), highlighting the capability of the proposed multi-variable calibration technique for improving the overall model performance. Further, the findings of this study highlighted the instrumental role of remote sensing data in representing the soil moisture dynamics of the study area and the importance of using multi-variable calibration to ensure robust hydrological simulations of river basins in dry climates.


Introduction
Water resource management in dry river basins is exceptionally challenging due to the competing stresses of prolonged droughts and sudden fluctuations in flow regimes, which can vary from low-flow conditions to devastating floods [1].Droughts, often extended in duration, lead to persistent water scarcity, compromising the availability of water for human consumption and the health of aquatic ecosystems [2].During these periods, water storage and supply systems are strained, reducing agricultural productivity and increasing competition among users, including cities, industries, and farmers.On the other hand, erratic rainfall results in sudden and extreme floods, which can overwhelm existing infrastructure, causing significant property damage, disrupting communities, and leading to potential loss of life.This variability necessitates better predictive models to anticipate and respond to changing conditions [3].
Hydrological models are inherently subject to uncertainties arising from model inputs, structure, and the quantification of initial conditions [4].As these uncertainties affect the reliability and accuracy of the model output, ultimately impacting sustainable water resource management, minimising them becomes a crucial consideration in hydrological modelling [5], which is usually achieved via optimising the model calibration.Although single-variable model calibration using streamflow data is the standard method for parameter optimisation, it often proves inadequate due to the complexity of various hydrological processes involved in water balance estimations [6].Accordingly, multi-variable or joint calibration techniques that combine soil moisture, evaporation data, and streamflow observations have proven effective [7][8][9].
In hydrological applications and simulations, soil moisture is a critical parameter due to its significant influence on various hydrological processes.Soil moisture affects rainfall partitioning into infiltration and runoff, directly impacting streamflow and groundwater recharge.Accurate soil moisture data, therefore, improve flood forecasting [10,11], drought predictions [12,13], wildfire mapping [14,15], irrigation scheduling [16,17], and the assessment of agricultural productivity [18,19].Conventional methods of obtaining soil moisture data rely on in situ point measurements collected from gauging stations.While continuous and accurate, these point measurements have inherent limitations, such as limited temporal coverage and the inability to represent regional soil moisture conditions within large river basins [20,21].
Advancements in remote sensing technology and satellite-based observations have significantly improved spatial and temporal coverage, offering extensive data even on subdaily time scales [21,22].For instance, the Soil Moisture and Ocean Salinity (SMOS) dataset from the European Space Agency and the Soil Moisture Active and Passive (SMAP) dataset from NASA provide global datasets that are widely used in hydrological studies.The SMOS dataset offers global coverage with spatial and temporal resolutions of 30~50 km and three days, respectively [23].The SMAP dataset has been available since March 2015, with spatial and temporal resolutions of 9 km and 3 h, respectively [22].Introducing such remotely sensed soil moisture observations to hydrologic model calibration has shown that joint calibration using soil moisture and streamflow observations enhances the robustness of model parameter projections, potentially leading to more accurate streamflow forecasts [21,24].
Among various hydrological models (viz., distributed (e.g., MIKE SHE [25] and TOP-MODEL [26]), semi-distributed (e.g., SWAT [27] and HEC-HMS [28]), and lumped (e.g., ABCD model [29] and HBV [30])), a lumped or conceptual model represents the aggregated behaviour of hydrological processes, allowing homogeneity and uniformity within the catchment.Despite their necessitated limitations in reflecting heterogeneous characteristics within the catchment, computational simplicity makes lumped hydrological models more effective in data-scarce regions [31].Frequently used lump hydrological models all over the world include the Tank Model [32] and the Hydrologiska Byråns Vattenbalansavdelning (HBV) Model [33], while the Hydrologic Engineering-Center Hydrologic Modelling System (HEC-HMS) with semi-distributed model capability is also often used for lump model applications.One commonly used lumped model is the bucket model, which partitions precipitation into several storage compartments.This approach reflects the land surface and sub-surface's capacity to store and transmit water based on a simple water balance concept [34].
The dry zone in Sri Lanka encompasses about 70% of the country and is characterised by low annual rainfall, typically below 1750 mm, concentrated during the monsoon season [35,36].The El Niño Southern Oscillation (ENSO) significantly influences the monsoon system in the Bay of Bengal, leading to inter-annual and intra-seasonal variability in monsoon rainfall [37,38].Consequently, persistent droughts in this region have severely impacted water resources and agricultural production [39].Additionally, these droughts have reduced food consumption, decreased investment capacity among farm households, and increased indirect costs such as healthcare expenses [40].Furthermore, climate projections for the 21st century indicate drastic changes in Sri Lanka, with increased rainfall in the wet zone and potentially severe droughts in the dry zone [41][42][43][44].In response to these changes, developing a methodology that better simulates basin hydrology is crucial for sustainable water resources management.Here, the Maduru Oya River Basin, one of the dry river basins in Sri Lanka, was selected to pilot the hydrological performance of a multi-variable model calibration method.In the initial phase of the study, the Tank Model parameters were calibrated and validated using observed streamflow data.The study then incorporates remote-sensing-based Root Zone Soil Moisture (RZSM) data from Soil Moisture Active Passive Level 4 (SMAP L4) and the observed streamflow to optimise the hydrological model parameters.The contribution of these two variables to the model parameter optimisation was estimated by introducing a weighted objective function, and the model performances were compared using the Kling-Gupta Efficiency.

Study Area
Maduru Oya is one of the major rivers in Sri Lanka, having a length of 135 km, and it flows through the intermediate and dry agricultural zones (Figure 1).The Maduru Oya catchment covers 1487 km 2 of land area, comprising a major tank (viz., Maduru Oya Reservoir) and six sub-watersheds (viz., upper, upper-middle, middle, and lower Maduru Oya, Kuda Oya, and Perillaveli Aru) [45].The Padiyathalawa sub-watershed, with a drainage area of 167 km 2 , is part of the upper Maduru Oya sub-catchment and is vital in supporting downstream agricultural activities.Scrublands and forests predominantly cover this sub-watershed.The primary soil types in the Padiyathalawa sub-watershed are reddish-brown earth and immature brown loams [46].

Hydro-Meteorological Data Collection and Preprocessing
The observed data required for model setup, calibration, and validation (viz., rainfall, maximum and minimum temperature, and streamflow) was obtained from the Department of Irrigation and the Department of Meteorology in Sri Lanka (Table 1).Considering the data availability and spatial distribution of rain gauges, precipitation data from Welipitiya Coconut Estate was selected for the study.Moreover, the missing data percentages The Padiyathalawa sub-watershed is located in the intermediate zone of Sri Lanka [47].It received an annual average precipitation of 2158 mm over the hydrological years between 1989 and 2020 [48].Despite this relatively sizeable amount of annual precipitation, the downstream of the Padiyathalawa River Basin (i.e., middle and downstream areas of the Maduru Oya basin) is in the dry zone.There is also significant seasonal variation in rainfall, with seven months (March to September) contributing less than 30% of the annual total rainfall, which includes two months (June and August) experiencing almost no rainfall (Figure A1).Significant seasonal variations in precipitation cause stress on numerous sectors, including agriculture and the ecological environment.The First Inter-Monsoon (March to April) and the South-West Monsoon (May to September) seasons bring light rains, contributing only ~9% and ~20% to the annual precipitation.In contrast, the Second Inter-Monsoon (October to November) and the North-East Monsoon (December to February) seasons account for ~29% and ~42% of the annual rainfall, respectively.The yearly average temperature of this sub-watershed is ~24

Hydro-Meteorological Data Collection and Preprocessing
The observed data required for model setup, calibration, and validation (viz., rainfall, maximum and minimum temperature, and streamflow) was obtained from the Department of Irrigation and the Department of Meteorology in Sri Lanka (Table 1).Considering the data availability and spatial distribution of rain gauges, precipitation data from Welipitiya Coconut Estate was selected for the study.Moreover, the missing data percentages for the surrounding precipitation gauging stations are less than 10%, thus making them suitable for the closest station patching method [49] to estimate missing climatic data.
Similarly, temperature records at Maduru Oya station, with less than 10% of missing data, were gap-filled using the closest station patching method.There are no missing data in the Padiyathalawa streamflow measurements over the study period considered.Monthly variations of each hydro-meteorological data (viz., rainfall, minimum and maximum temperature, and streamflow) are shown in Appendix A Figures A1-A4.
Evaporation is one of the inputs required by the Tank Model.Since there are no in-situ measurements available, Hargreaves' reference evapotranspiration model (1985) was used to calculate the reference evapotranspiration (ET 0 ) at the Maduru Oya temperature gauge station (Figure 1).The limitations of a simple evapotranspiration model can be addressed by locally calibrating it using ET values obtained from a pan evaporation method [50].Accordingly, the ET 0 values were calibrated with the in situ long-term average monthly pan evaporation (pan coefficient of 0.8) recorded at the Padiyathalawa gauging station [48].The computation of Reference Evapotranspiration (ET 0 ) is presented in Appendix B, and the calibration results are presented in Table A1.
Remote-sensing-based soil moisture data (freely available at the National Snow and Ice Data Centre (NSIDC) [22]) and observed streamflow were used for the model calibration and validation.As the SMAP L4 retrievals are obtainable at the 3 h frequency for a particular location, it provides eight RZSM values for a given day.Preliminary data checking found no outliers in the dataset during the simulation period.Therefore, due to the minimal daily variations in these data, they were averaged using the simple arithmetic mean method to obtain daily values.The temporally aggregated data was then spatially aggregated using the Thiessen polygon method.It was noticed that the spatially and temporally aggregated RZSM data from SMAP comprise systematic biases.Even though the remote sensing data capture temporal variation, the Tank Model simulated soil moisture data show differences in its magnitude.Therefore, the spatially and temporally aggregated RZSM data was rescaled using the Cumulative Distribution Function (CDF) Matching method, and its outcome is shown in Figure A5.

Hydrological Model Set-Up
In this study, we selected the Tank Model, a lumped, non-linear, process-oriented Conceptual Rainfall-Runoff model developed by Sugawara [32].Here, the series storage type is considered the best model among different Tank Models available (viz., exponential, parallel exponential, overflow, storage, and series storage) [32].This model consists of several tanks vertically connected in series with bottom and side outlets to replicate the distinct watershed hydrological process (Figure 2).The series storage type Tank Model simulates the runoff as a function of the precipitation, evaporation, and antecedent soil moisture (or water storage in the soil), resulting in a non-linear deterministic model characteristic [52].This model structure analyses the daily streamflow of an outlet consuming daily precipitation and evaporation inputs.The non-linear structure of the model accounts for initial abstraction [53].It distributes the rainfall and assigns time lags to each component automatically.Further, it represents a reasonable physical representation corresponding to the groundwater structure.Despite its capability to describe the non-linear nature of surface runoff, it is challenging to determine the model parameters due to this non-linear behaviour [53].Moreover, depending on the modelling objective and catchment characteristics, the series-storage type Tank Model has variations in the number of tanks and hydrological processes primarily considered.Among the numerous studies employing the Tank Model, Gunawardhana and Kazama [54] used it for a low-flow analysis under climate change conditions in the Tagliamento River in Italy.Yokoo et al. [55] regionalised the Tank Model parameters by relating them to catchment characteristics such as land use, topography, geology, and soil type.Goodazri et al. [56] investigated the effect of catchment scale on the Tank Model parameters, concluding that the model can be effectively applied even in larger river basins.
Many successful studies in which the Tank Model was used in different regions can be found in the literature (e.g., [53,[57][58][59]).The Tank Model has been successfully employed in catchments as small as 15 km 2 (Goodazri et al. [56]) and as large as 1935 km 2 (Gunawardhana and Kazama [54]).Following Yokoo et al. [55] and Musiake and Wijesekera [60], the initial model parameters were estimated using catchment characteristics such as land use, topography, geology, and soil type.Model parameters were then further calibrated using observed river flow data to match the specific site conditions in the study area.When focusing on high flows, the parameters directly related to the surface runoff, such as A 11 , Z 11 , and A 12 , are identified as the most significant.On the other hand, low-flow characteristics are governed by the sub-baseflow and baseflow parameters, which respond to the river flow with a time lag.The variables of the hydrological model are determined by the site-specific characteristics.For example, in snow-fed river basins, a snow melt component governed by air temperature is added to the first tank [54].Considering long-term simulations, evapotranspiration was included in the water balance computations of this study.As this study focuses on a dry watershed with seasonal variations in precipitatio the Tank Model developed for humid regions is applied [53].The typical hydrologi model parameters are directly measured physical or calibrated process parameters.Ho ever, the Tank Model is designed for the calibrated process parameters [62].The Ta Model parameters can be calibrated either manually (using a trial-and-error process) [6 or automatically by the use of computer-aided complex mathematical formulations [6 Although automatic calibration improves the efficiency, stability, and consistency of t model performance [63,65], its major drawback is the determination of initial paramet that govern the computer algorithm to achieve the optimal parameters [52].General determining the realistic initial value becomes challenging because it requires experien and user knowledge of the model and the study area [52].However, the effectiveness the Tank Model depends on the degree of calibration [45], which can be improved by u ing automatic calibration in conjunction with manual calibration [57,65].Initial model p rameters used in this study (Table 2) were assigned based on the Tank Model paramet optimised for the Mahaweli River Basin in Sri Lanka [60].
A sensitivity analysis was conducted to identify the parameter sensitivity to the pe and low flow simulation, in which the adopted initial values varied between −75% a 75% to evaluate their contributions to the observed flow between October 2010 and Se tember 2015 (Figure 3).A one-year warm-up period was considered before two datas were considered during the study to initialise the water storage in each tank under thr cycles, assuming all tanks were empty during the initial round.The model parameters required in the Tank Model can be clustered into two groups: (1) outlet coefficient (five parameters on the side walls, viz., A 11 , A 12 , A 2 , A 3 , and A 4 and three parameters on the tank bottoms, viz., B 1 , B 2 , and B 3 ) and ( 2) soil water storage coefficients (four runoff hole depths, viz., Z 11 , Z 12 , Z 2 , and Z 3 ) [52,57,59,61].The output from the side outlet of the top tank (i.e., Tank 1) is surface runoff, and the output from the second tank is intermediate runoff.The outputs from the third and fourth tanks are considered to account for low-flow situations via the sub-base flow and baseflow, respectively (Figure 2).Mathematical interpretations of the Tank Model are illustrated in the following Equations ( 1)-( 4), where x is the tank count, starting from the top; n is the count of days since the simulation starts (days); ∆t is the time interval (days); A x is the runoff coefficient of the x th tank (1/day); B x is the infiltration coefficient of the x th tank (1/day); H x is the water depth in the x th tank (mm); Z x is the runoff hole depth of the x th tank (mm); Q x is the contribution to the streamflow from the x th tank (mm/day); I x is the infiltration from the x th tank (mm/day); P is precipitation (mm/day); E is evapotranspiration (mm/day), and Q T is the total streamflow (mm/day).As this study focuses on a dry watershed with seasonal variations in precipitation, the Tank Model developed for humid regions is applied [53].The typical hydrological model parameters are directly measured physical or calibrated process parameters.However, the Tank Model is designed for the calibrated process parameters [62].The Tank Model parameters can be calibrated either manually (using a trial-and-error process) [63] or automatically by the use of computer-aided complex mathematical formulations [64].Although automatic calibration improves the efficiency, stability, and consistency of the model performance [63,65], its major drawback is the determination of initial parameters that govern the computer algorithm to achieve the optimal parameters [52].Generally, determining the realistic initial value becomes challenging because it requires experience and user knowledge of the model and the study area [52].However, the effectiveness of the Tank Model depends on the degree of calibration [45], which can be improved by using automatic calibration in conjunction with manual calibration [57,65].Initial model parameters used in this study (Table 2) were assigned based on the Tank Model parameters optimised for the Mahaweli River Basin in Sri Lanka [60].Note: Initial parameter values were taken from Musiake and Wijesekera [60]; the minimum and maximum values of the model parameters were taken from Sugawara et al. [53].
A sensitivity analysis was conducted to identify the parameter sensitivity to the peak and low flow simulation, in which the adopted initial values varied between −75% and 75% to evaluate their contributions to the observed flow between October 2010 and September 2015 (Figure 3).A one-year warm-up period was considered before two datasets were considered during the study to initialise the water storage in each tank under three cycles, assuming all tanks were empty during the initial round.

Objective Function for Model Optimisation
The performance of hydrological model simulations can be evaluated by numerous statistical indices, including the Nash-Sutcliffe efficiency (NSE), Nash and Sutcliffe, Kling-Gupta efficiency (KGE) [66], and Root Mean Square Error (RMSE).The Nash-Sutcliffe Efficiency (NSE) tends to be overly sensitive to high-flow values.Moreover, its use of the observed mean as a baseline can result in overestimating model performance for highly seasonal variables, such as river flow in dry basins [67].One of the major limitations of the RMSE is that it is more influenced by a few large errors in the simulation, which may not provide a balanced view of model performance [68].Both NSE and RMSE provide a single aggregated measure of error, which can mask the underlying reasons for poor model performance.
The KGE was explicitly designed to balance the sensitivity to bias and variability.It decomposes the error into three components: correlation (r), variability (α), and bias (β), allowing an understanding of different aspects of model performance, indicating temporal pattern, mean difference, and dispersion difference between observed and simulated data, respectively.Therefore, the KGE ensures a more balanced evaluation of the model performance across different flow conditions [69].
The Kling-Gupta Efficiency is defined as the Euclidean Distance (ED) to the point measured from the ideal point (5).The KGE values can vary between −∞ and 1.0, with the latter indicating a perfect model accuracy [66]. 2  (5) Here, r is the linear correlation coefficient between the observed and simulated data (6), α is the relative variability between observed and simulated data (7), and β denotes the bias (8); where µ S and µ O are the mean values of simulated and observed data, respectively, σ S and σ S are the standard deviations of simulated and observed data, respectively.

Model Calibration and Optimisation
The model calibration procedure involves two steps, single-and multi-variable, undertaken between August 2011 and September 2012.The calibrated model was validated with data between September 2015 and April 2016.

Single-Variable Calibration
The single-objective calibration approach was employed to optimise the model with observed daily streamflow data utilising KGE (KGE Q ) for the calibration period.The model was calibrated based on the initial parameters developed for the upper Peradeniya catchment (Table 2) [60].The results were validated using the observed data between September 2015 and April 2016.

Multi-Variable Calibration
A weighted objective function (9) was used under 11 calibration schemes, where the weighted factor (ϑ) varies between 0 and 1.0, with 0.1 increments.This calibration was undertaken between September 2015 and April 2016, constrained by the availability of calibration variables (viz., streamflow and remote-sensing-based soil moisture).
Microsoft Excel's Solver includes one linear and two non-linear solver functions (Generalized Reduced Gradient (GRG) and evolutionary solver), and the GRG solver has been widely used in hydrological studies.Given the non-linear nature of the hydrological processes involved in basin water balance, this study employs the GRG method to estimate the parameters of the Tank Model, focused on optimising the weighted KGE function (9).All the twelve model parameters (Table A2) were set as changing variables, constrained by their respective upper limits.For multi-variable optimisation, the GRG Non-linear method [70] was selected while assuming the optimisation problem to be smooth and non-linear, which is an assumption involved in other studies (e.g., [71][72][73]).
Here, the Multi Start option was enabled for a population size of 20, which commences the calibration with unique starting points per each population class.Since the random seed was set to zero, the program selects each starting point randomly, preserving the diversity of solutions.The chosen population size yielded a satisfactory convergence (tolerance margin of 0.0001) through single-objective optimisation, swiftly converging the multi-variable optimisation to an adequate level.Furthermore, two commonly used tests (viz., t-test [74] and Wilcoxon signed-rank test [75]) were employed to identify the statistical significance of the model performance improvement in two calibration schemes.

Results
The content of this section is structured into three sub-sections.The first sub-section presents the findings of the parameter sensitivity analysis.The second and third subsections present performances under single-and multi-variable calibrations.

Parameter Sensitivity
The results of the sensitivity analysis are presented in Figure 3. Here, the model parameters can be categorised based on their sensitivity levels towards peak and low flows (highly sensitive, moderately sensitive, and less sensitive).The highest flow rate during the simulation period was used as the reference value for peak flow analysis, while the 25th percentile of the flow rate was used for low-flow analysis.A higher percentile was selected for the low-flow sensitivity analysis due to the large number of zero flow rates (no flows) occurring during the dry period.Parameters of the first tank (A 11 , A 12 , Z 11 , Z 12 , and B 1 ), which are related to infiltration and surface, sub-surface, and intermediate flow, as well as the runoff coefficient of the second tank (A 2 ), were found to have higher sensitivity towards high flows.In contrast, parameters in other tanks (Z 2 , B 2 , Z 3 , A 3 , B 3 , and A 4 ) were found to be trivial for producing high flow rates.These less-sensitive parameters represent the intermediate flow, infiltration, base flow, and sub-base flow of different tanks.Baseflow and sub-baseflow parameters were found to be more sensitive for low flows.For example, increases in parameters Z 2 and B 3 positively affect low flows, while decreases in parameters A 11 , A 2 , and A 3 also positively affect low flows.This phenomenon can be explained by the fact that as surface runoff decreases due to lower corresponding parameters, infiltration and percolation increase, enhancing low flows.

Model Performances in Single-Variable Calibration
The single-variable calibration using observed streamflow data provides satisfactory streamflow simulations in both calibration (Figure 4) and validation (Figure 5), with respective KGE Q values of 0.87 and 0.85.However, the hydrographs (Figure 4A) and the Flow Duration Curves (FDCs) (Figure 4B) in calibration show marginal underestimations in the peak-flow and low-flow simulations.

Model Performances in Single-Variable Calibration
The single-variable calibration using observed streamflow data provides satisfactory streamflow simulations in both calibration (Figure 4) and validation (Figure 5), with respective  values of 0.87 and 0.85.However, the hydrographs (Figure 4A) and the Flow Duration Curves (FDCs) (Figure 4B) in calibration show marginal underestimations in the peak-flow and low-flow simulations.

Model Performances in Multi-Variable Calibration
The multi-variable calibration performances in terms of KGE are presented in Figure 6.The two extreme calibration schemes, viz.,  = 0 and  = 1.0 represent the single variable optimisations with streamflow and root zone soil moisture, respectively.Except for the  = 1 calibration scheme (i.e., single variable optimisation with soil moisture), model performances in streamflow and soil moisture show satisfactory results under all other

Model Performances in Multi-Variable Calibration
The multi-variable calibration performances in terms of KGE are presented in Figure 6.The two extreme calibration schemes, viz., ϑ = 0 and ϑ = 1.0 represent the single variable optimisations with streamflow and root zone soil moisture, respectively.Except for the ϑ = 1 calibration scheme (i.e., single variable optimisation with soil moisture), model performances in streamflow and soil moisture show satisfactory results under all other schemes tested in this study.The calibration scheme of ϑ = 1 shows that the single objective optimisation with soil moisture fails to produce accurate streamflow simulations.Figure 6 clearly shows that the weighted KGE (KGE W , ( 9)) increases from ϑ = 0.2 scheme onwards, causing slight increases in KGE SM and slight decreases in KGE Q with a drastic decrease at ϑ = 1 scheme.Therefore, it is evident that, from a multi-variable perspective, higher weightage on root zone soil moisture can enhance the overall model performance in terms of KGE W . Figure 7 illustrates the percentage change in  resulting from multi-variable calibration compared to single-variable across all calibration schemes.Among the accuracyincreased calibration schemes (ϑ values of 0.1, 0.8, and 0.9), the ϑ = 0.1 scheme yields the most accurate streamflow simulation, with a  of 0.87 and 0.89 in calibration and validation, respectively.However, despite offering a better simulation with  of 0.90, it exhibits lower accuracy in soil moisture simulation than other schemes (Figure 6).Regarding the calibration schemes, the best RZSM simulation is observed with the ϑ = 0.9 scheme.Nonetheless, since this study's primary focus is on accurate streamflow simulation in the catchment, the calibration scheme of ϑ = 0.1 is deemed more suitable as it provides the highest accuracy in streamflow simulation, even during the validation period with a  of 0.89. Figure 7 illustrates the percentage change in KGE Q resulting from multi-variable calibration compared to single-variable across all calibration schemes.Among the accuracyincreased calibration schemes (ϑ values of 0.1, 0.8, and 0.9), the ϑ = 0.1 scheme yields the most accurate streamflow simulation, with a KGE Q of 0.87 and 0.89 in calibration and validation, respectively.However, despite offering a better simulation with KGE SM of 0.90, it exhibits lower accuracy in soil moisture simulation than other schemes (Figure 6).Regarding the calibration schemes, the best RZSM simulation is observed with the ϑ = 0.9 scheme.Nonetheless, since this study's primary focus is on accurate streamflow simulation in the catchment, the calibration scheme of ϑ = 0.1 is deemed more suitable as it provides the highest accuracy in streamflow simulation, even during the validation period with a KGE Q of 0.89.The simulated streamflow and FDCs for the better-performing (i.e., accuracy-increased) calibration schemes with weighted factor (ϑ) of 0.1, 0.8, and 0.9 are illustrated in Figure 8.It shows that the calibration scheme of ϑ = 0.1 underestimates the peak and low flows.However, compared to the single-variable calibration (i.e., ϑ = 0), the ϑ = 0.1 calibration scheme provides better peak and lower flow estimations.The other two calibration schemes (i.e., ϑ = 0.8 and ϑ = 0.9) considerably underestimate the low flows (Figure 8).
The streamflow simulations and FDCs of the better-performing schemes during the validation period (between August 2011 and September 2012) are shown in Figure 9.It shows that, even for low flow conditions, the calibration scheme ϑ = 0.1 provides better estimations than the single-variable calibration (i.e., ϑ = 0).As identified in the calibration, the calibration schemes ϑ = 0.8 and ϑ = 0.9 cannot simulate low flow conditions satisfactorily.
The soil moisture variation under accuracy-increased calibration schemes of ϑ = 0.1, 0.8, and 0.9 are shown in Figure 10.As discussed earlier, the calibration scheme ϑ = 0.1 shows the minimum agreement between the observed and simulated soil moisture.The calibration scheme with a 0.9 weighted factor (ϑ) shows better agreement in observed and simulated soil moisture content.The simulated streamflow and FDCs for the better-performing (i.e., accuracy-increased) calibration schemes with weighted factor (ϑ) of 0.1, 0.8, and 0.9 are illustrated in Figure 8.It shows that the calibration scheme of ϑ = 0.1 underestimates the peak and low flows.However, compared to the single-variable calibration (i.e., ϑ = 0), the ϑ = 0.1 calibration scheme provides better peak and lower flow estimations.The other two calibration schemes (i.e., ϑ = 0.8 and ϑ = 0.9) considerably underestimate the low flows (Figure 8).
The streamflow simulations and FDCs of the better-performing schemes during the validation period (between August 2011 and September 2012) are shown in Figure 9.It shows that, even for low flow conditions, the calibration scheme ϑ = 0.1 provides better estimations than the single-variable calibration (i.e., ϑ = 0).As identified in the calibration, the calibration schemes ϑ = 0.8 and ϑ = 0.9 cannot simulate low flow conditions satisfactorily.
The soil moisture variation under accuracy-increased calibration schemes of ϑ = 0.1, 0.8, and 0.9 are shown in Figure 10.As discussed earlier, the calibration scheme ϑ = 0.1 shows the minimum agreement between the observed and simulated soil moisture.The calibration scheme with a 0.9 weighted factor (ϑ) shows better agreement in observed and simulated soil moisture content.

Discussion
Our results show increased model performance (i.e., better KGE values) under multivariable calibration than a single variable.As an example, calibration schemes with weighted factors (ϑ) of 0.1, 0.8, and 0.9 show more accurate streamflow simulations than the single-variable calibration, even in the validation period.These better performances can be attributed to the behaviour of different calibration parameters.For example, A11 and A12 coefficients, which enhance the surface and sub-surface flows, are higher for the ϑ = 0.1 than ϑ =0 (viz., single-variable calibration with streamflow), resulting in more accurate simulations of the peak flows.Also, higher B1, B2, and B3 values enhance the water percolation up to the groundwater table, and greater A3 and A4 coefficients increase the base-flow contribution.Even though overall model performances ( ) are higher under ϑ = 0.8 and 0.9 calibration schemes, those schemes significantly underestimate the low flows due to considerably low A3 and A4 coefficients (Table A2) (compared to the same in calibration under ϑ = 0.1).Previous studies (e.g., Xiong and Zeng [21]) have concluded that using multiple variables enhances the reliability of optimised parameters, even if the overall model performance is less than that of single-variable optimisation.The KGE used in this study assesses different aspects of model performance, including correlation, variability, and bias.We evaluated the increases in model performance using the paired twosample t-test and Wilcoxon signed-rank test, in which 8.5% and 24.7% statistical improvements are observed with the respective tests.
The literature shows that the multi-variable calibration using in situ streamflow and remote sensing-based soil moisture data can marginally reduce the streamflow simulation during the calibration compared with single-variable calibration with streamflow [21,76].Further, Xiong and Zeng [21] utilised SMAP Level 3 surface soil moisture retrievals converted to the Soil Wetness Index of the root zone, potentially resulting in a less accurate representation of catchment soil-water conditions compared to the direct RZSM retrievals from SMAP L4 used in this study.The conclusions drawn by Gunasekara et al. [76] are

Discussion
Our results show increased model performance (i.e., better KGE values) under multivariable calibration than a single variable.As an example, calibration schemes with weighted factors (ϑ) of 0.1, 0.8, and 0.9 show more accurate streamflow simulations than the single-variable calibration, even in the validation period.These better performances can be attributed to the behaviour of different calibration parameters.For example, A 11 and A 12 coefficients, which enhance the surface and sub-surface flows, are higher for the ϑ = 0.1 than ϑ =0 (viz., single-variable calibration with streamflow), resulting in more accurate simulations of the peak flows.Also, higher B 1 , B 2 , and B 3 values enhance the water percolation up to the groundwater table, and greater A 3 and A 4 coefficients increase the base-flow contribution.Even though overall model performances (KGE Q ) are higher under ϑ = 0.8 and 0.9 calibration schemes, those schemes significantly underestimate the low flows due to considerably low A 3 and A 4 coefficients (Table A2) (compared to the same in calibration under ϑ = 0.1).Previous studies (e.g., Xiong and Zeng [21]) have concluded that using multiple variables enhances the reliability of optimised parameters, even if the overall model performance is less than that of single-variable optimisation.The KGE used in this study assesses different aspects of model performance, including correlation, variability, and bias.We evaluated the increases in model performance using the paired two-sample t-test and Wilcoxon signed-rank test, in which 8.5% and 24.7% statistical improvements are observed with the respective tests.
The literature shows that the multi-variable calibration using in situ streamflow and remote sensing-based soil moisture data can marginally reduce the streamflow simulation during the calibration compared with single-variable calibration with streamflow [21,76].Further, Xiong and Zeng [21] utilised SMAP Level 3 surface soil moisture retrievals converted to the Soil Wetness Index of the root zone, potentially resulting in a less accurate representation of catchment soil-water conditions compared to the direct RZSM retrievals from SMAP L4 used in this study.The conclusions drawn by Gunasekara et al. [76] are controversial, even though their research utilised the lumped hydrological model (ABCD Model, ref. [29]) and SMAP L4 RZSM data in a daily time resolution similar to this study.However, Gunasekara et al.'s study [76] differs from ours in that it features a significant water body within the study area, which poses challenges to the applicability of a lumped hydrological model.
Simply put, the Tank Model used in this study successfully conceptualised the hydrology of the study area.Eini et al. [77] also used streamflow and Climate Change Initiative Soil Moisture dataset in the multi-variable calibration of the Soil and Water Assessment Tool+ (SWAT+, ref. [78]) in a sizeable transboundary basin (viz., Odra River Basin) in Central Europe.They concluded that the multi-objective calibration could substantially improve the accuracy of the simulations with an average KGE of 0.67 compared to an average KGE of 0.60 in single-variable calibration.A study by Loizu et al. [24] used Advanced SCATterometer Surface Soil Moisture data to improve the streamflow simulations in Mediterranean catchments through different data assimilation techniques.Therefore, careful evaluation of potential remote sensing soil moisture products is essential in hydrological model calibration to improve the performance of streamflow and other water balance components simulations in varied regions.
Besides soil moisture data, other remote-sensing-based variables are used in multivariable calibration to enhance the hydrological model performances, particularly in datascarce regions.For example, Sirisena et al. [8] used remote sensing-based evapotranspiration data from Global Land Evaporation: the Amsterdam Model (GLEAM ET) along with streamflow in the Chindwin River Basin in Myanmar, where in situ measurements are sparse.They also suggested that multi-variable calibration offers good performance in streamflow while maintaining reasonable performances when simulating the catchment's evapotranspiration.Furthermore, Rientjes et al. [79] indicated that the water balance of the Karkheh River Basin, Iran is best reproduced by streamflow and satellite-based evapotranspiration data used during the hydrological model calibration.Therefore, it is evident that the potential use of different satellite-based products as alternative sources in hydrological model calibration may return dividends over traditional calibration, particularly in data-poor regions with limited streamflow records.

Conclusions
This study focused on assessing the significance of multi-variable calibration in hydrological simulations in a data-scarce dry river basin by utilising in situ streamflow and remote sensing-based soil moisture data (viz., SMAP L4) in single-and multi-variable calibration schemes using a lumped hydrological mode (viz., the Tank Model).The proposed method was piloted at the Padiyathalawa sub-basin of Madura Oya Basin, Sri Lanka.For the multi-objective calibration, a series of weighted factors (ϑ), between 0 and 1.0 (at 0.1 intervals), were used, and the performances were assessed using the Kling-Gupta Efficiency (KGE).This study is one of the few studies conducted in the region to investigate the possibility of improving hydrological model performances using Soil Moisture Active Passive Level 4 (SMAP L4) under different calibration schemes (i.e., single-and multivariable).Furthermore, the lumped Tank Model can successfully represent the catchment dynamics, particularly in small catchments such as the Padiyathalawa sub-watershed.
With single-variable calibration (i.e., streamflow only), both calibration and validation showed higher model efficiency criteria (KGE > 0.83), even though it slightly underestimates the peak flows while considerably underestimating the low flows.Three multi-variable calibration schemes with 0.1, 0.8, and 0.9 weighted factors (ϑ) showed better streamflow simulations (KGE Q > 0.85 for both calibration and validation) than the single-variable calibration.Therefore, incorporating root zone soil moisture into the weighted multi-objective function improves the overall performance of the Tank Model.However, the calibration scheme of 0.1 weighted factor (ϑ) is more capable of simulating the streamflow dynamics, particularly the low flows, than the other two better-performing calibration schemes (ϑ of 0.8 and 0.9).Therefore, the calibration scheme with a 0.1 weighted factor (ϑ) is more suitable for conceptualising the catchment dynamics at the Padiyathalawa sub-watershed, underscoring the effect of the minimal contribution from root zone soil moisture translating into significantly improved streamflow simulations in terms of accuracy and reliability.
It is important to note that this study's calibration and validation periods were intentionally kept short due to significant inconsistencies in the hydro-meteorological data.As the calibration and validation periods within this limited dataset may not accurately represent the long-term catchment flow dynamics, it is essential to validate these results with long-term data better to assess the robustness of the proposed modelling technique.Also, exploring alternative multi-variable optimisation techniques, such as genetic algorithms, is recommended.Furthermore, it is worth exploring the suitability of other remotesensing-based soil moisture products and hydrological variables like evapotranspiration in enhancing the hydrological model's performance within the study area.

Figure 2 .
Figure 2. Conceptual illustration of the Tank Model setup.

Figure 2 .
Figure 2. Conceptual illustration of the Tank Model setup.

Figure 3 .
Figure 3. Sensitivity analysis of the Tank Model parameters for (A) peak flow and (B) first quartile (i.e., 25th percentile) flow.Note: A11 -runoff coefficient to estimate surface flow, Z11 -runoff hole depth to estimate surface flow, A12 -runoff coefficient to estimate sub-surface flow, Z12 -runoff hole depth to estimate sub-surface flow, B1 -infiltration coefficient of the first tank, A2 -runoff coefficient to estimate intermediate flow, Z2 -runoff hole depth to estimate intermediate flow, B2 -infiltration coefficient of the second tank, A3 -runoff coefficient to estimate sub-base flow, Z3 -runoff hole depth to estimate sub-base flow, B3 -infiltration coefficient of the third tank, A4 runoff coefficient to estimate base flow.

Figure 3 .
Figure 3. Sensitivity analysis of the Tank Model parameters for (A) peak flow and (B) first quartile (i.e., 25th percentile) flow.Note: A 11 -runoff coefficient to estimate surface flow, Z 11 -runoff hole depth to estimate surface flow, A 12 -runoff coefficient to estimate sub-surface flow, Z 12 -runoff hole depth to estimate sub-surface flow, B 1 -infiltration coefficient of the first tank, A 2 -runoff coefficient to estimate intermediate flow, Z 2 -runoff hole depth to estimate intermediate flow, B 2 -infiltration coefficient of the second tank, A 3 -runoff coefficient to estimate sub-base flow, Z 3 -runoff hole depth to estimate sub-base flow, B 3 -infiltration coefficient of the third tank, A 4 runoff coefficient to estimate base flow.

Figure 4 .
Figure 4. Comparison of observed and simulated hydrographs (A) and flow-duration curves (B) at the basin outlet for the calibration period (between August 2011 and September 2012).

Figure 4 . 25 Figure 5 .
Figure 4. Comparison of observed and simulated hydrographs (A) and flow-duration curves (B) at the basin outlet for the calibration period (between August 2011 and September 2012).Hydrology 2024, 11, x FOR PEER REVIEW 12 of 25

Figure 5 .
Figure 5.Comparison of observed and simulated hydrographs (A) and flow-duration curves (B) at the basin outlet for the validation period (between September 2015 and April 2016).

Hydrology 2024 , 25 Figure 7 .
Figure 7. Percentage change in Kling-Gupta Efficiency for streamflow ( ) when the optimal parameters from multi-variable calibration are used in the simulation, compared to the simulations with parameters optimised from single-variable calibration (i.e., weighted factor (ϑ) of 0).

Figure 7 .
Figure 7. Percentage change in Kling-Gupta Efficiency for streamflow (KGE Q ) when the optimal from multi-variable calibration are used in the simulation, compared to the simulations with parameters optimised from single-variable calibration (i.e., weighted factor (ϑ) of 0).

Funding:
This research received no external funding.

Figure A1 .
Figure A1.Monthly precipitation at Welipitiya Coconut Estate gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the horizontal line shows the median (i.e., 50th percentile) value of the monthly data sets.Whiskers are extended to 1.5 times inter-quartile range to the top and bottom of the boxes.

Figure A2 .
Figure A2.Average minimum temperature measured at Maduru Oya gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the horizontal line shows the median (i.e., 50th percentile) value of the monthly data sets.Whiskers are extended to 1.5 times inter-quartile range to the top and bottom of the boxes.

Figure A3 .
Figure A3.Average maximum temperature measured at Maduru Oya gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the

Figure A2 .
Figure A2.Average minimum temperature measured at Maduru Oya gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the horizontal line shows the median (i.e., 50th percentile) value of the monthly data sets.Whiskers are extended to 1.5 times inter-quartile range to the top and bottom of the boxes.

Figure A1 .
Figure A1.Monthly precipitation at Welipitiya Coconut Estate gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the horizontal line shows the median (i.e., 50th percentile) value of the monthly data sets.Whiskers are extended to 1.5 times inter-quartile range to the top and bottom of the boxes.

Figure A2 .
Figure A2.Average minimum temperature measured at Maduru Oya gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the horizontal line shows the median (i.e., 50th percentile) value of the monthly data sets.Whiskers are extended to 1.5 times inter-quartile range to the top and bottom of the boxes.

Figure A3 .
Figure A3.Average maximum temperature measured at Maduru Oya gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the

Figure A3 .
Figure A3.Average maximum temperature measured at Maduru Oya gauging station between October 2010 and September 2016.The boxes are limited to the 25th and 75th percentiles, and the horizontal line shows the median (i.e., 50th percentile) value of the monthly data sets.Whiskers are extended to 1.5 times inter-quartile range to the top and bottom of the boxes.

Table 1 .
In situ and remote sensing data used for the analysis.

Table 2 .
Initial model parameters and their ranges.
: * Kling-Gupta Efficiency (KGE) over the validation period (August 2011-September 2012) and the rest are the KGE values over the calibration period (September 2015-April 2016).KGE Q , KGE SM , and KGE W are the streamflow, soil moisture, and weighted Kling-Gupta Efficiency, respectively.B 1 -infiltration coefficient of the first tank, A 11 -runoff coefficient to estimate surface flow, A 12 -runoff coefficient to estimate sub-surface flow, Z 11 -runoff hole depth to estimate surface flow, Z 12 -runoff hole depth to estimate sub-surface flow, B 2 -infiltration coefficient of the second tank, A 2 -runoff coefficient to estimate intermediate flow, Z 2 -runoff hole depth to estimate intermediate flow, B 3 -infiltration coefficient of the third tank, A 3 -runoff coefficient to estimate sub-base flow, Z 3 -runoff hole depth to estimate sub-base flow, A 4 -runoff coefficient to estimate base flow. Note