Dam Water Level Prediction Using Vector AutoRegression, Random Forest Regression and MLP-ANN Models Based on Land-Use and Climate Factors

: To predict the variability of dam water levels, parametric Multivariate Linear Regression (MLR), stochastic Vector AutoRegressive (VAR), Random Forest Regression (RFR) and Multilayer Perceptron (MLP) Artiﬁcial Neural Network (ANN) models were compared based on the inﬂuences of climate factors (rainfall and temperature), climate indices (DSLP, Aridity Index (AI), SOI and Niño 3.4) and land-use land-cover (LULC) as the predictor variables. For the case study of the Gaborone dam and the Bokaa dam in the semi-arid Botswana, from 2001 to 2019, the prediction results showed that the linear MLR is not robust for predicting the complex non-linear variabilities of the dam water levels with the predictor variables. The stochastic VAR detected the relationship between LULC and the dam water levels with R 2 > 0.95; however, it was unable to sufﬁciently capture the inﬂuence of climate factors on the dam water levels. RFR and MLP-ANN showed signiﬁcant correlations between the dam water levels and the climate factors and climate indices, with a higher R 2 value between 0.890 and 0.926, for the Gaborone dam, compared to 0.704–0.865 for the Bokaa dam. Using LULC for dam water predictions, RFR performed better than MLP-ANN, with higher accuracy results for the Bokaa dam. Based on the climate factors and climate indices, MLP-ANN provided the best prediction results for the dam water levels for both dams. To improve the prediction results, a VAR-ANN hybrid model was found to be more suitable for integrating LULC and the climate conditions and in predicting the variability of the linear and non-linear time-series components of the dam water levels for both dams.


Introduction
Despite freshwater scarcity being a global problem, the solutions must be locally formulated in order to understand the connections between water supply and demand, and to adequately respond to the local water shortages. Water shortages are being exacerbated by human activities, as manifested by population growth and the impacts of land-use changes driven by urbanization, agricultural activities, industrialization, and economic development. The cumulative effects of the intensification of land-use activities and climate change continue to pose uncertainty on the availability of water resources, with the effect of intensified manipulation of the surface and groundwater hydrological regimes [1].

Study Area
The study area is located within Botswana's Limpopo River Basin (BLRB). The larger Limpopo River Basin is a transboundary basin, covering an area of approximately 416,300 km 2 , and straddles four southern African countries: South Africa (45%), Botswana (19%), Mozambique (21%), and Zimbabwe (15%). The basin is home to more than 18 million people and Botswana has the highest percentage (61%) of its population living in the basin. As shown in Figure 1, the semi-arid Botswana relies on the following small-tomedium-sized dams, which are located within the BLRB: Gaborone (141. 4  of this work is on the derivation of a hybrid model capable of combining stochastic and machine learning models for the accurate prediction of dam water levels in the two dams by integrating the LULC and the climate conditions within the dam catchments.

Study Area
The study area is located within Botswana's Limpopo River Basin (BLRB). The larger Limpopo River Basin is a transboundary basin, covering an area of approximately 416,300 km 2 , and straddles four southern African countries: South Africa (45%), Botswana (19%), Mozambique (21%), and Zimbabwe (15%). The basin is home to more than 18 million people and Botswana has the highest percentage (61%) of its population living in the basin. As shown in Figure 1, the semi-arid Botswana relies on the following small-to-mediumsized dams, which are located within the BLRB: Gaborone (141. 4   With the general scarcity of freshwater in the arid and semi-arid regions, water management problems tend to worsen, especially during extreme hydrological events, such as drought. For this reason, and to optimally manage the dam operations, continuous and accurate reservoir management schemes-including predictions of the variabilities of the dam water capacities and the determinations of the influences of a natural climatic phenomenon and anthropogenic activities on the water resources-is essential. In most regions, predicting and forecasting dam water capacities is still challenging for water resource operators and managers. This is attributed to the fact that, despite reservoir water levels being directly regulated by the inflows and outflow releases, there are several uncertainties in the dam water level determinant variables, such as the temporal dynamics of climatic factors, e.g., rainfall and temperature, and dam operations and management regimes, which are complex to model. With the general scarcity of freshwater in the arid and semi-arid regions, water management problems tend to worsen, especially during extreme hydrological events, such as drought. For this reason, and to optimally manage the dam operations, continuous and accurate reservoir management schemes-including predictions of the variabilities of the dam water capacities and the determinations of the influences of a natural climatic phenomenon and anthropogenic activities on the water resources-is essential. In most regions, predicting and forecasting dam water capacities is still challenging for water resource operators and managers. This is attributed to the fact that, despite reservoir water levels being directly regulated by the inflows and outflow releases, there are several uncertainties in the dam water level determinant variables, such as the temporal dynamics of climatic factors, e.g., rainfall and temperature, and dam operations and management regimes, which are complex to model. atmospheric correction algorithm and the Landsat rescaling coefficients, the multitemporal Landsat images were corrected to generate the surface reflectance imagery. The LULC classification was carried out using Breiman's random forest algorithm [28] and was implemented within the Google Earth Engine, as detailed in [29]. To improve the classification accuracy, the mean, variance, and contrast gray-level cooccurrence matrix (GLCM) texture features were found to be most significant and were included in the classification scheme. The LULC classification accuracy metrics results are presented in Table 1, and the LULC area coverages are summarized in Table 2 for the Bokaa and Gaborone dam catchments. From the results, the Bokaa dam catchment occupies an area of approximately 3610 km 2 and the Gaborone catchment is approximately 4344 km 2 . From the classification error matrix, the overall accuracy (OA) is determined from the ratio of the correctly classified pixels to the total training sample. Further, the respective class User's Accuracy (UA) is determined by the ratio of the correct positive predictions, while the Producer's Accuracy (PA) is the ratio of the correctly detected positives. For each year, the average OA, UA and PA are presented in Table 1. The results in Table 1 show that for both dam catchments, the LULC classification accuracies, as measured using the PA, UA and OA metrics, were higher than 80%, and the corresponding Kappa Index ranged between 0.75 and 0.87. The accuracy measures demonstrate that the LULC was derived with a high degree of accuracy for both dam catchments. In Table 2, it is observed that for both catchments, the built-up areas are increasing  exponentially, while the vegetation and bare soil-covered areas increased and decreased, interchangeably, either due to activities in croplands or due to climate influences. Tree cover within the catchments is also observed to be increasing in coverage, while shrubland has decreased in extent over the years.

1.
Rainfall and Temperature Monthly rainfall data from the Gaborone gauge station was also used for both the Bokaa dam and Gaborone dam catchments due to their geographical proximity, climatic similarities, and given that there is no gauge station within the Bokaa dam catchment. Figure 2 shows the observed rainfall patterns within the Bokaa and Gaborone dam catchments, and Figure 3 shows the minimum, average, and maximum temperature variabilities within the catchments. Over the 19 years of study, it is observed that the mean temperature is increasing while the amount of rainfall received in the two catchments is decreasing. Monthly rainfall data from the Gaborone gauge station was also used for both the Bokaa dam and Gaborone dam catchments due to their geographical proximity, climatic similarities, and given that there is no gauge station within the Bokaa dam catchment. Figure 2 shows the observed rainfall patterns within the Bokaa and Gaborone dam catchments, and Figure 3 shows the minimum, average, and maximum temperature variabilities within the catchments. Over the 19 years of study, it is observed that the mean temperature is increasing while the amount of rainfall received in the two catchments is decreasing.

Climate Indices
The climate indices considered were those that have teleconnections with particular rainfall over southern Africa, that is, DSLP, SOI, and Niño 3.4. The average March-June pressures at Darwin have proven to have high positive sea level pressure (SLP) anomalies and teleconnections to droughts over southern Africa [30]. The SOI standardized sea-level pressure difference between Papeete and Darwin is also related to rainfall over the subregion. In addition to the three climate indices, the aridity index (AI) was derived using station rainfall and temperature data, as in Equation (1): 12 i P

Climate Indices
The climate indices considered were those that have teleconnections with particular rainfall over southern Africa, that is, DSLP, SOI, and Niño 3.4. The average March-June pressures at Darwin have proven to have high positive sea level pressure (SLP) anomalies Sustainability 2022, 14, 14934 7 of 31 and teleconnections to droughts over southern Africa [30]. The SOI standardized sea-level pressure difference between Papeete and Darwin is also related to rainfall over the subregion. In addition to the three climate indices, the aridity index (AI) was derived using station rainfall and temperature data, as in Equation (1): where P i = the monthly total precipitation (mm) and T i = mean near-surface temperature ( • C).

Dam Reservoir Water Levels
The mean monthly dam water levels were used as the indicators for water availability in surface water storages, from 2001 to 2019, for the two dams. Figure 4a shows the variability of the dam water levels with rainfall, with the Bokaa dam exhibiting a marginally higher degree of correlation with rainfall than the Gaborone dam. The scatterplot regressions in Figure 4b depicts very low correlations between the measured dam water levels and rainfall in the two catchments.

Climate Indices
The climate indices considered were those that have teleconnections with particular rainfall over southern Africa, that is, DSLP, SOI, and Niño 3.4. The average March-June pressures at Darwin have proven to have high positive sea level pressure (SLP) anomalies and teleconnections to droughts over southern Africa [30]. The SOI standardized sea-level pressure difference between Papeete and Darwin is also related to rainfall over the subregion. In addition to the three climate indices, the aridity index (AI) was derived using station rainfall and temperature data, as in Equation (1): where i P = the monthly total precipitation (mm) and i T = mean near-surface temperature (°C).

Dam Reservoir Water Levels
The mean monthly dam water levels were used as the indicators for water availability in surface water storages, from 2001 to 2019, for the two dams. Figure 4a shows the variability of the dam water levels with rainfall, with the Bokaa dam exhibiting a marginally higher degree of correlation with rainfall than the Gaborone dam. The scatterplot regressions in Figure 4b depicts very low correlations between the measured dam water levels and rainfall in the two catchments.

Data Statistics and Correlational Analysis
The summary of the mean monthly statistical descriptions of the study datasets, from 2001 to 2019, for the two dams is presented in Table 3.
In terms of the correlations presented in Figure 5, the Bokaa dam exhibits the highest water levels, but inverse correlations with tree cover, at −0.349, followed by maximum temperature, bare soil, grassland, average temperature and aridity index, respectively, at −0.243, 0.216, 0.175, 0.161, and 0.161. The Bokaa dam water level correlations were particularly worse with Niño 3.4 and DSLP, at −0.035 and −0.047, respectively. In general, the water levels in the Bokaa dam have positive but low correlations with LULC classes and low negative correlations with the climate factors. Comparatively, the Gaborone dam had higher correlations with the predictor variables ( Figure 5). The highest correlations for the Gaborone dam water levels with the predictor variables were for grassland, water bodies, and shrubland, at 0.815, 0.761 and −0.730, respectively. The lowest correlations were with built up, rainfall, and aridity index, at 0.013, 0.029 and 0.034, respectively. The Gaborone dam displays positive and higher correlations with dam surface area and grassland; however, lower and negative correlations with climate factors and indices. MLR was utilized as a baseline for competition evaluation [24]. Linear regression models are simple models that have linear and non-linear parameters for predictions. For small sample sizes, the parametric multilinear regression (MLR) models are able to establish the relationships between the predictor variables and the dependent variable using least squares fitting. In this study, the dam water levels depend on climate factors, climate indicators, and LULC. The general MLR model is expressed as in Equation (2).

Multivariate Linear Regression (MLR)
MLR was utilized as a baseline for competition evaluation [24]. Linear regression models are simple models that have linear and non-linear parameters for predictions. For small sample sizes, the parametric multilinear regression (MLR) models are able to establish the relationships between the predictor variables and the dependent variable using least squares fitting. In this study, the dam water levels depend on climate factors, climate indicators, and LULC. The general MLR model is expressed as in Equation (2).

Vector AutoRegressive Model
VAR is a stochastic linear prediction model that predicts the current time variable value, based on its previous time value, and takes into consideration other predictor variables. Through dynamic analysis, VAR detects the changes to a particular variable, affects changes to other variables, the lags of those variables and the changes in the variables' lags. VAR thus extends the univariate autoregression to the multiple time-series regression, with the lagged values of all series as regressors. For example, the VAR model of two

Vector AutoRegressive Model
VAR is a stochastic linear prediction model that predicts the current time variable value, based on its previous time value, and takes into consideration other predictor variables. Through dynamic analysis, VAR detects the changes to a particular variable, affects changes to other variables, the lags of those variables and the changes in the variables' lags. VAR thus extends the univariate autoregression to the multiple time-series regression, with the lagged values of all series as regressors. For example, the VAR model of two variables X t and Y t (k = 2) with the lag order p is defined as in Equations (3) and (4). The β and γ can be estimated using the ordinary least squares method.
The lag-order for the VAR(p) model is determined using the lag-length selection criteria, and the VAR(p) models are fitted with orders p = 0, 1, . . . p max and the p-value, which minimizes some model selection criteria, is chosen. The parameter lag selection criteria in this study are the Akaike's Information Criterion (AIC p ), Schwarz Bayesian Information Criterion (BIC p ), Hannan-Quinn Criterion (HQC p ), and Final Prediction Error (FPE p ). The traditional unrestricted VAR is unsuitable for non-stationary data with seasonality and, therefore, this study imposed a priori differencing on the input datasets for stationarity. The implemented VAR model for dam water level time series prediction was developed with the following steps:

1.
Testing for stationarity of the individual predictor variables using the augmented Dickey-Fuller (ADF) test.

2.
Determining the lag for the VAR(p) model using lag-length selection. VAR(p) models are fitted with orders p = 0, 1, . . . p max , and the p value resulting in minimal model selection criteria is chosen based on the parameter selection criteria above. In this study, the lag orders are determined for the specific predictor variables.

3.
Establishment of an optimal VAR model with appropriate lags for each parameter. For multivariate time series, the VAR model is constructed such that each variable, at a time point, exhibits as a linear function of the recent lag of itself and other variables. The generalized VAR(p) = VAR(1) form for the n = 15 predictor variables can be expressed as in Equation (5). Equation (5) is solved using ordinary least squares, and c represents the intercepts; A is the regression coefficient matrix, and e is the error in prediction at time t.
Residual autocorrelation assessment for goodness-of-fit. For the time series data, the autocorrelation of the residuals between the observed and the model-fitted values is used to determine the goodness-of-fit of the model. Accuracy assessment metrices, including R 2 , RMSE, MAE and MAPE are used. 5.
VAR system stability test assessment with the autoregressive (AR) roots graph. The VAR stability determines how well the model represents the time series over the sampling window. This is evaluated using the roots of the characteristic polynomial of the coefficient matrix A in Equation (5). If the roots are less than 0, the VAR model is considered stable.

Random Forest Regression
RFR is an ensemble learning regression model based on a decision tree algorithm [28]. The RFR principle entails randomly generating different unpruned CART decision trees, in which the decrease in Gini impurity is regarded as the splitting criterion. As a bootstrap resampling and bagging approach, the bootstrap samples from the training set data are fitted with an unpruned decision tree for each bootstrap sample. At the decision tree nodes, variable selection is made on small random subsets of the predictor variables and the best split from the predictors used to split the node. The trees in the forest are averaged or voted to generate output probabilities and a final model that generates a robust model. In this study the construction of the RFR through the following steps:

1.
From the original data, nTree bootstrap samples are drawn.

2.
For each bootstrap dataset, a tree is grown, and for each tree-node mTry variables are randomly selected for splitting.

3.
The aggregated information from the nTree trees is used for new data prediction, in this case voting for regression.

4.
Out-of-bag (OOB) error rate are computed using the test dataset not in the bootstrap sample.
RFR hyperparameters were tuned to determine the optimal lag-order, epochs, number of trees (n_estimators) and max_depth for predicting the dam water levels.

Multilayer Perceptron (MLP) Neural Network
MLP-ANN is one of the most popular Neural Network models with input, hidden, and output layers. The advantage of MLP-ANN is that even with a single hidden layer and arbitrary bounded and smooth activation function, the network can approximate a continuous non-linear system. The adopted network in this study was trained on the Levenberg-Marquardt backpropagation with a gradient scheme for weighting adjustment to minimize the predicted and observed data errors. The MLP-ANN model was implemented following the structure and detailed steps in [31].

Performance Evaluation Metrics
Four statistical measures were used to evaluate the prediction efficiency of the models, RMSE, R 2 , MAE, and MAPE. The metrices are respectively represented in Equations (6) where h o i is the observed dam water level and h s i is the simulated or predicted dam water level. RMSE, MAE, and h are measured in % of dam water level.

Data Normalization
The input datasets were standardized within the range [0.1-0.9]. The [0.1-0.9] normalization, using the minimum-maximum boundary, was used to standardize the original data, as expressed in Equation (10). The standardization minimizes biases as all the input data receive the same attention. To evaluate the relative importance of the predictor variables, backward sensitivity analysis is adopted, where the significance of each input variable is determined by stepwise variable replacement and the measure of the MAE deviation.

Parameter Lag Order Determination for VAR Model
The optimal lag orders for the Gaborone and Bokaa dams were determined based on the AIC, BIC, and HQIC measures. From the summary results in Table 4, rainfall (Set-4) had the lowest AIC, BIC, and HQIC information criteria for the Gaborone dam, respectively corresponding to −9.493, −5.990, and −8.070 (Table 4). Set-7, comprising rainfall and temperature, was the second lowest, followed by Set-8, consisting of all the temperatures, and the highest measure was detected from Set-3, comprising all the parameters. For the Gaborone dam, temperature and rainfall had the highest lag orders, at 43 and 40, respectively. From the results in Table 4 for the Bokaa dam, the rainfall factor (Set-4) gave the lowest AIC, BIC, and HQIC, at −9.061, −7.732 and −8.523, respectively, and the highest lag order of 20. This is followed by Set-7, combining rainfall and temperatures, with a lag order of 12. For the Bokaa dam, the respective optimal lag orders varied between 7-20, with temperature having the least lag order compared to the Gaborone dam. In Table 4, the FPE values are not included since their magnitudes were all negligible. The VAR training results show that the contributions of rainfall and temperatures were insignificant for both dams, with R 2 of less than 35%. The combination of the two climate factors in Set-7 only improved the training results for the Bokaa dam water levels, but did not influence the water levels in the Gaborone dam. Both dams responded well with LULC and the four regional climate indices, with R 2 of between 76% and 95% (Table 4).

Training for RF Regression
To determine the optimal RFR tuning hyperparameters, the data sets were trained with 70% of the data. The training results, based on lag order and max_depth, n_estimators, are presented in Figure 6 for the Bokaa dam, and the corresponding results for the four best predictors variables are presented in Table 5. For the Gaborone dam water level simulations and predictions, the results for the RFR model tuning parameters are also presented in Figure 6, with the best predictor variables statistics presented in Table 5. The RFR hyperparameter tuning results show that the water level prediction in the Bokaa dam required significantly higher lag orders than the Gaborone dam but relatively shallower depth and fewer n_estimators or number of RFR trees ( Table 5). The RFR training results for the best datasets depict R 2 > 0.82, with the exception of the Gaborone dam, where the climate indices yield R 2 = 0.563.  The RFR hyperparameter tuning results show that the water level prediction in the Bokaa dam required significantly higher lag orders than the Gaborone dam but relatively shallower depth and fewer n_estimators or number of RFR trees ( Table 5). The RFR training results for the best datasets depict R 2 > 0.82, with the exception of the Gaborone dam, where the climate indices yield R 2 = 0.563.

Training of MLP-ANN Model
The training of the MLP-ANN for predicting the water levels in the two dams was based on the lag order, the network number of hidden layers, epochs, and batch sizes. The tuning results for the dams are illustrated in Figure 7 and the summary statistics for the best four predictor variables are presented in Table 6.
For MLP training, low lag orders, between 1-4, are required to train the ANN, with the hidden layers varying from 2-4 ( Table 6). The Bokaa dam required higher epochs, with relatively lower batch sizes, to train the model compared to the Gaborone dam, with the exception of the data set comprising min-avg-max temperatures for the Gaborone dam. The difference between the MLP and RFR hyperparameter tuning is that MLP-ANN detected the direct impact of rainfall (Set-4) on the Bokaa dam water level variability, while RFR only detected it indirectly, in combination with temperature (Set-7). For the Gaborone dam, RFR detected the direct impact of climate indices (Set-8), however, this was only captured indirectly for the Bokaa dam using RFR with Set-1. The RFR and MLP-ANN results indicate that the temporal variability of the dam water levels within the two catchments is influenced by the climate indices and climate factors. The impact of LULC is not directly related to the water levels but may contribute to the determination of demand and dam operation regimes. The best predictor variables in Table 6 show high training output with R 2 > 0.83.
The observed variable responses in the hyperparametric tuning for water levels in both dams, using RFR and MLP-ANN, respectively, shown in Figures 6 and 7, are attributed to the systematic one-parameter-at-a-time tuning approach. For both models, the input of the combination of the predetermined optimal hyperparameters in the determination of the final hyperparameter automatically minimizes the model errors yielding the best fit results, as observed in the final tuning response curves in For MLP training, low lag orders, between 1-4, are required to train the ANN, with the hidden layers varying from 2-4 ( Table 6). The Bokaa dam required higher epochs, with relatively lower batch sizes, to train the model compared to the Gaborone dam, with the exception of the data set comprising min-avg-max temperatures for the Gaborone dam. The difference between the MLP and RFR hyperparameter tuning is that MLP-ANN detected the direct impact of rainfall (Set-4) on the Bokaa dam water level variability, while RFR only detected it indirectly, in combination with temperature (Set-7). For the Gaborone dam, RFR detected the direct impact of climate indices (Set-8), however, this was only captured indirectly for the Bokaa dam using RFR with Set-1. The RFR and MLP-ANN results indicate that the temporal variability of the dam water levels within the two catchments is influenced by the climate indices and climate factors. The impact of LULC is not directly related to the water levels but may contribute to the determination of demand and dam operation regimes. The best predictor variables in Table 6 show high training output with R 2 > 0.83.

Dam Water Level Prediction Results
This section presents the standardized dam water prediction results for comparison between the two dams. The RMSE, MAE, and MAPE are calculated on the inverse of Equation (9) of the standardized datasets.

Prediction of Dam Water Levels Using MLR
The results for the prediction of the dam water levels using MLR are presented in Table 7, and shows that for both dams, Set-3, comprising of all variables, Set-5 (LULC), and Set-6 (LULC, Rainfall, Min and Max Temperatures) were the best predictors. For the Bokaa dam, the highest R 2 was 0.583, from Set-3 and Set-6, while the same sets yielded R 2 = 0.841 for the Gaborone dam, and LULC (Set-5) had R 2 of 0.785 for the Gaborone dam, compared to 0.489 for the Bokaa dam. The rest of the predictor variables predicted the time-series variability of the dam water levels at less than 50% accuracy in terms of R 2 . Since the same regression fitting equation was used for training and testing the time-series dam water levels, the MLR results were found to be similar, with very low prediction accuracy. Using the same fit for the entire 19-year data gave better results, as presented in Figure 8, and demonstrated the fact that more robust model(s), at both training and testing phases, are required in the prediction of dam water levels. The dam water level predictions for the Bokaa dam were based on the predetermined optimal training results for each dataset, shown in Table 4. The prediction results show that only Sets-1, -3, -5, -6, and -8 presented the highest convergence for the Bokaa dam (Table 8)  Despite the good predictions using LULC for the Gaborone dam, which impacts on Set-3 and Set-6, the graphical plots in Figure 8 and the large RMSE, MASE, and MAPE show that the linear MLR is not suitable for simulating and predicting the complex, seasonal, and non-linear trends exhibited by the water levels in both dams. As such, the MLR results confirm the hypothesis that more robust regression models are necessary for predicting water levels in the dams.

Bokaa Dam Water Level Prediction Using VAR
The dam water level predictions for the Bokaa dam were based on the predetermined optimal training results for each dataset, shown in Table 4. The prediction results show that only Sets-1, -3, -5, -6, and -8 presented the highest convergence for the Bokaa dam (Table 8). Set-5, comprising the LULC classes, gave the highest R 2 , at 0.998. The second highest (R 2 = 0.975) predictor variable is (Set-6), followed by Set-1 (R 2 = 0.959), Set-3 (R 2 = 0.928) and Set-8 (R 2 = 0.916). In terms of climate indices and climate factors, Set-1 (R 2 = 0.959; RMSE = 3.3%; MAE = 2.7%; MAPE = 14.3%) and Set-8 (R 2 = 0.995; RMSE = 2.7%; MAE = 2.2%; MAPE = 36.9%) gave the best results. Without the climate indices, the longterm predictions of dam water levels using temperatures (Set-2), rainfall (Set-4), and their combination shows low prediction results. The rainfall and temperature sets registered the highest MAPE errors, of more than 50%. The good performance of the LULC is attributed to the interpolation within the five years, which results in minimal variability within the input data and, therefore, low data variability and high accuracy.  The best results for the water level predictions in the Bokaa dam are presented in Figure 9. The prediction results and the graphical fits show that, despite having the highest performance accuracy, the predictor factors combined LULC (Set-3 and Set-6) are not the best predictor variables. This is particularly due to the inability of the model to capture the dam water levels at the beginning of the prediction using the LULC as the predictor factor. These differences are captured within the dotted boxes in Figure 9, depicting a lack of expected trends and patterns. From the graphical and statistical analysis, the best predictor variables for the Bokaa dam water levels are Set-1 and Set-8, where Set-1 was influenced by both climate factors and climate indices.

Gaborone Dam Water Level Prediction Using VAR
Using the VAR model, the prediction of the Gaborone dam water levels is detected to be significant using the four climate indices (Set-8), as shown in Table 8 and Figure 10 (R 2 = 0.929; RMSE = 0.7%; MAE = 0.6%; MAPE = 8%). The rainfall and temperature climate factors performed marginally in predicting the Gaborone dam water levels, with R 2 of less than 0.3 and MAPE above 20%, while their combination in Set-7 yielded higher accuracy prediction accuracy results. Similarly, high prediction results were obtained using the integration of the climate indices with rainfall and temperature in Set-1. The results for the climate-based predictors are presented in Figure 10 for the Gaborone dam, with Set-3 including all parameters. By visually assessing the trends of the predictions within the dotted boxes in Figure 10, it is empirically observed that climate indices gave the best results. However, the results show that in the absence of climate factors, LULC can be used to predict the water levels in the dams with good accuracy (R 2 > 0.990; RMSE < 0.7%; MAE < 0.3%; MAPE < 3.5%). est performance accuracy, the predictor factors combined LULC (Set-3 and Set-6) are not the best predictor variables. This is particularly due to the inability of the model to capture the dam water levels at the beginning of the prediction using the LULC as the predictor factor. These differences are captured within the dotted boxes in Figure 9, depicting a lack of expected trends and patterns. From the graphical and statistical analysis, the best predictor variables for the Bokaa dam water levels are Set-1 and Set-8, where Set-1 was influenced by both climate factors and climate indices. Figure 9. Water level prediction for Bokaa dam using VAR model. Reprinted with permission from ref. [27]. Copyright 2022 Society of Photo-Optical Instrumentation Engineers.

Gaborone Dam Water Level Prediction Using VAR
Using the VAR model, the prediction of the Gaborone dam water levels is detected to be significant using the four climate indices (Set-8), as shown in Table 8 and Figure 10 (R 2 = 0.929; RMSE = 0.7%; MAE = 0.6%; MAPE = 8%). The rainfall and temperature climate factors performed marginally in predicting the Gaborone dam water levels, with R 2 of less than 0.3 and MAPE above 20%, while their combination in Set-7 yielded higher accuracy prediction accuracy results. Similarly, high prediction results were obtained using the integration of the climate indices with rainfall and temperature in Set-1. The results for the climate-based predictors are presented in Figure 10 for the Gaborone dam, with Set-3 including all parameters. By visually assessing the trends of the predictions within the dotted boxes in Figure 10, it is empirically observed that climate indices gave the best results. However, the results show that in the absence of climate factors, LULC can be used to predict the water levels in the dams with good accuracy (R 2 > 0.990; RMSE < 0.7%; MAE < 0.3%; MAPE < 3.5%).

RFR Prediction of Bokaa Dam Water Levels
The RFR prediction results for the Bokaa dam show that all the datasets are suitable for predicting the water levels, with R 2 > 0.8. LULC and RFR presented the least prediction accuracy with R 2 = 0.807 and the best four predictors were Set-2 of all the temperatures, followed by Set-3, Set-7, and Set-1, with R 2 of 0.836, 0.829, 0.824, and 0.820, respectively.

RFR Prediction of Bokaa Dam Water Levels
The RFR prediction results for the Bokaa dam show that all the datasets are suitable for predicting the water levels, with R 2 > 0.8. LULC and RFR presented the least prediction accuracy with R 2 = 0.807 and the best four predictors were Set-2 of all the temperatures, followed by Set-3, Set-7, and Set-1, with R 2 of 0.836, 0.829, 0.824, and 0.820, respectively. The corresponding RMSE varied between 11.3-12.5%, with an MAE average of approximately 7% and MAPE of approximately 13%. Figure 11 presents the predictions for the four best predictor variable sets. The results from Set-2 and Set-7 comprise temperatures and rainfall and depict that RFR captured the relationship between the dam water levels and the climate factors (rainfall and temperature). The analysis of the prediction trends confirms Set-2 and Set-7 as the most suitable for predicting the dam water levels, as illustrated within the dotted boxes, where the predictor variables are able to capture the temporal trends of the measured dam water levels.

RFR Prediction of Gaborone Dam Water Levels
Using the optimal RFR hyperparameters for predicting the water levels in the Gaborone dam, Table 9 shows Sets-2, -4, -7, and -8 presenting the best results, with R 2 values of 0.918, 0.819, 0.898, 0.897 and 0.890, respectively. The datasets comprise temperature, rainfall, their combination, and the climate indices, respectively. The RMSE is observed to be lower for the Gaborone dam than the Bokaa dam, ranging between 9.7% and 11.4%, while the MAE averages were at 6.5% of dam water levels and MAPE is higher, at between 23% and 38%. The LULC-based prediction results show that, despite the positive correlation of more than 65%, with the dam water levels, LULC does not capture the temporal seasonality and variability of the dam water levels ( Figure 12). The results in Table 9 and Figure 12 depict that RFR is able to predict the water levels in the Gaborone dam using the climate factors, with the temperatures (Set-2) being the best climate factor, followed by rainfall (Set-4). The combination of temperature and rainfall marginally reduces the influence of the predictive ability of temperatures by nearly 10%, to R 2 of 0.898. The climate indices (Set-8) display a significant impact on the water levels in the Gaborone dam, with R 2 = 0.890. The dotted box regions in Figure 12 show the inability of RFR to accurately Figure 11. Observed and RFR predicted water levels for Bokaa dam. Reprinted with permission from ref. [27]. Copyright 2022 Society of Photo-Optical Instrumentation Engineers.

RFR Prediction of Gaborone Dam Water Levels
Using the optimal RFR hyperparameters for predicting the water levels in the Gaborone dam, Table 9 shows Sets-2, -4, -7, and -8 presenting the best results, with R 2 values of 0.918, 0.819, 0.898, 0.897 and 0.890, respectively. The datasets comprise temperature, rainfall, their combination, and the climate indices, respectively. The RMSE is observed to be lower for the Gaborone dam than the Bokaa dam, ranging between 9.7% and 11.4%, while the MAE averages were at 6.5% of dam water levels and MAPE is higher, at between 23% and 38%. The LULC-based prediction results show that, despite the positive correlation of more than 65%, with the dam water levels, LULC does not capture the temporal seasonality and variability of the dam water levels ( Figure 12). The results in Table 9 and Figure 12 depict that RFR is able to predict the water levels in the Gaborone dam using the climate factors, with the temperatures (Set-2) being the best climate factor, followed by rainfall (Set-4). The combination of temperature and rainfall marginally reduces the influence of the predictive ability of temperatures by nearly 10%, to R 2 of 0.898. The climate indices (Set-8) display a significant impact on the water levels in the Gaborone dam, with R 2 = 0.890. The dotted box regions in Figure 12 show the inability of RFR to accurately predict the temporal trends in the Gaborone dam water levels.

Bokaa Dam Water Level Prediction Using MLP-ANN
With the rectifier linear unit activation function, Adam optimizer, and a learning rate of 0.0003, the results for predicting water levels in the Bokaa dam are presented in Table  10. The local temperature is linked to the dam water levels with the highest R 2 of 0.865, and the lowest RMSE = 10.9% and MAE = 6.5%. The combination of temperature and rainfall (Set-7) is second, with R 2 of 0.850, followed by rainfall (R 2 = 0.829). Climate indices (Set-8) also influenced the dam water levels with R 2 of 0.805 and the least MAPE = 13.2%. LULC had the least influence on the dam water levels, with MAPE of 27.7%, and its combination with the other parameters in Set-6 further reduced the accuracy, with MAPE = 56.6% and R 2 = 0.449.

Bokaa Dam Water Level Prediction Using MLP-ANN
With the rectifier linear unit activation function, Adam optimizer, and a learning rate of 0.0003, the results for predicting water levels in the Bokaa dam are presented in Table 10. The local temperature is linked to the dam water levels with the highest R 2 of 0.865, and the lowest RMSE = 10.9% and MAE = 6.5%. The combination of temperature and rainfall (Set-7) is second, with R 2 of 0.850, followed by rainfall (R 2 = 0.829). Climate indices (Set-8) also influenced the dam water levels with R 2 of 0.805 and the least MAPE = 13.2%. LULC had the least influence on the dam water levels, with MAPE of 27.7%, and its combination with the other parameters in Set-6 further reduced the accuracy, with MAPE = 56.6% and R 2 = 0.449. The performances for the best predictor variables within the box window time regions in Figure 13 show that local temperature (Set-2) and rainfall (Set-4) exhibit similar and best prediction trends with MAEs of approximately 6.5% and MAPE of 25%. The performances for the best predictor variables within the box window time regions in Figure 13 show that local temperature (Set-2) and rainfall (Set-4) exhibit similar and best prediction trends with MAEs of approximately 6.5% and MAPE of 25%.

Gaborone Dam Water Level Prediction Using MLP-ANN
For the Gaborone dam, all the predictor datasets with LULC (Sets-3-5-6) did not converge to predict the dam water levels (Table 10). This further confirms the observations in MLR and RFR, where LULC recorded low correlations with dam water levels. The best performing sets in predicting dam water levels for the Gaborone dam were Set-4, rainfall (0.926), performing equally with Set-2 (0.925), then Set-7 (0.920), and Set-1 (0.917). The results show a positive response of the dam water levels to rainfall, temperature, and to the climate indices with an average low RMSE of less than 10%, R 2 > 0.91, and the least MAE, > 5% on average. The dotted boxes in Figure 14 show the differences in the dam water predictions for the Gaborone dam. In comparison to Set-7, Sets-1, -2, and -4 present good initial estimations of the dam water level. The MLP-ANN results improved the ability of RFR to detect near-linear trends, with Set-2 (temperature) presenting the best em-

Gaborone Dam Water Level Prediction Using MLP-ANN
For the Gaborone dam, all the predictor datasets with LULC (Sets-3-5-6) did not converge to predict the dam water levels (Table 10). This further confirms the observations in MLR and RFR, where LULC recorded low correlations with dam water levels. The best performing sets in predicting dam water levels for the Gaborone dam were Set-4, rainfall (0.926), performing equally with Set-2 (0.925), then Set-7 (0.920), and Set-1 (0.917). The results show a positive response of the dam water levels to rainfall, temperature, and to the climate indices with an average low RMSE of less than 10%, R 2 > 0.91, and the least MAE, >5% on average. The dotted boxes in Figure 14 show the differences in the dam water predictions for the Gaborone dam. In comparison to Set-7, Sets-1, -2, and -4 present good initial estimations of the dam water level. The MLP-ANN results improved the ability of RFR to detect near-linear trends, with Set-2 (temperature) presenting the best empirical and statistical predictions (Figure 14).

Relative Importance of the Predictor Variables
For the Bokaa dam, Figure 15 presents the relative importance of each variable in the predictor groups and compares all the factors. Comparing the variables, the tree-cover and shrubland exhibited the highest correlation with dam water levels (slightly more than 50% influence), followed by the max temperature and Niño 3.4. The least contributions are from bare soil, built-up and grassland, with the significance of rainfall and aridity index being negligible. The significance of the predictor variables indicates that within the Bokaa catchment, the degree of vegetation index and the regional temperature have higher correlations with the Bokaa dam water capacity. For the Gaborone dam ( Figure  15), grassland and water bodies exhibit the highest significance, followed by cropland and bare soil, with the rest of the parameters contributing less than 2% each. The aridity index and rainfall are observed to have the least contributions toward predicting the Gaborone dam water levels. While grassland has negligible contributions to dam water levels in the Bokaa dam, it has the highest significance for water capacity in the Gaborone dam, accounting for nearly 48% significance. Similar to the Bokaa dam, the significance of vegetation health is observed to have higher correlations with the dam water levels in the Gaborone dam.

Relative Importance of the Predictor Variables
For the Bokaa dam, Figure 15 presents the relative importance of each variable in the predictor groups and compares all the factors. Comparing the variables, the tree-cover and shrubland exhibited the highest correlation with dam water levels (slightly more than 50% influence), followed by the max temperature and Niño 3.4. The least contributions are from bare soil, built-up and grassland, with the significance of rainfall and aridity index being negligible. The significance of the predictor variables indicates that within the Bokaa catchment, the degree of vegetation index and the regional temperature have higher correlations with the Bokaa dam water capacity. For the Gaborone dam (Figure 15), grassland and water bodies exhibit the highest significance, followed by cropland and bare soil, with the rest of the parameters contributing less than 2% each. The aridity index and rainfall are observed to have the least contributions toward predicting the Gaborone dam water levels. While grassland has negligible contributions to dam water levels in the Bokaa dam, it has the highest significance for water capacity in the Gaborone dam, accounting for nearly 48% significance. Similar to the Bokaa dam, the significance of vegetation health is observed to have higher correlations with the dam water levels in the Gaborone dam. Investigating the predictor data groups for the Bokaa dam, in terms of the catchme LULC, tree-cover has the most influence in predicting the dam water levels, accountin for more than 50%; built-up, bare-soil, and grassland have the least contribution, with th significance of each at less than 1%. The climate factors and maximum temperature hav the highest contributions, at 34%, and rainfall at 22% for the Bokaa catchment. Among th climate indicators, Niño 3.4 has the highest contribution in predicting the dam water le els in the Bokaa dam, at 28%. For the Gaborone dam, the existence of water bodies an grassland is most important for predicting the dam water levels in the Gaborone dam with up to 32%. The climate factors exhibit competing significance, ranging between 21 25%, with minimum temperature and rainfall as the most significant climate factors. F the climate indices, Niño 3.4 has the highest significance, at 42%, with AI and DSLP bein the least, with a nearly equal relative importance of 17%.
The relative importance measures, shown in Figure 16, depict the sensitives of th predictor variables. The results show that for both dams, LULC forms part of the mo significant predictor variables; therefore, the more accurate catchment LULC, in terms high temporal resolution and actual classification accuracy, is important in predictin dam water levels for both dams. The parametric sensitivities in Figures 15 and 16 al Investigating the predictor data groups for the Bokaa dam, in terms of the catchment LULC, tree-cover has the most influence in predicting the dam water levels, accounting for more than 50%; built-up, bare-soil, and grassland have the least contribution, with the significance of each at less than 1%. The climate factors and maximum temperature have the highest contributions, at 34%, and rainfall at 22% for the Bokaa catchment. Among the climate indicators, Niño 3.4 has the highest contribution in predicting the dam water levels in the Bokaa dam, at 28%. For the Gaborone dam, the existence of water bodies and grassland is most important for predicting the dam water levels in the Gaborone dam, with up to 32%. The climate factors exhibit competing significance, ranging between 21-25%, with minimum temperature and rainfall as the most significant climate factors. For the climate indices, Niño 3.4 has the highest significance, at 42%, with AI and DSLP being the least, with a nearly equal relative importance of 17%.
The relative importance measures, shown in Figure 16, depict the sensitives of the predictor variables. The results show that for both dams, LULC forms part of the most significant predictor variables; therefore, the more accurate catchment LULC, in terms of high temporal resolution and actual classification accuracy, is important in predicting dam water levels for both dams. The parametric sensitivities in Figures 15 and 16 also imply that the prediction model should be able to capture the influences of both the high and low significant variables. imply that the prediction model should be able to capture the influences of both the high and low significant variables.

Discussions
The present study compares the performance of the stochastic VAR and the machine learning RFR and MLP-ANN models. The performances of each prediction horizon are compared using the average MAE, RMSE, and MAPE estimates and the R 2 statistics as a goodness-of-fit of the models. The metrics are considered to adequately measure the prediction accuracy and depict how well the model generalizes the unseen or test data. To determine the best predictor variables and to gauge the sensitivity of the models to the inputs, different exogenous input combinations were explored, and the results were compared using the above statistical indicators.

Impact of LULC on Water Level Predictions
The current study reveals the significance of LULC in predicting dam water levels as detected by the tested models. The assumption in the five-year time epoch used in the LULC temporal resolution is that there are insignificant changes in the natural land-covers such as water bodies, grasslands, shrublands, forests, bare soils, and land-use such as croplands. However, significant changes are expected in urban built-up, although at a slow spatial and temporal rate. Only the stochastic VAR detected the correlation and variability between the dam water levels and LULC, and predicted the dam water levels with LULC as the best predictor variable, with the highest accuracy of greater than 99%. The prediction results using MLR, RFR, and MLP-ANN showed that the LULC pattern, as interpolated over the 20-year period, may not be suitable for predicting the dam water

Discussions
The present study compares the performance of the stochastic VAR and the machine learning RFR and MLP-ANN models. The performances of each prediction horizon are compared using the average MAE, RMSE, and MAPE estimates and the R 2 statistics as a goodness-of-fit of the models. The metrics are considered to adequately measure the prediction accuracy and depict how well the model generalizes the unseen or test data. To determine the best predictor variables and to gauge the sensitivity of the models to the inputs, different exogenous input combinations were explored, and the results were compared using the above statistical indicators.

Impact of LULC on Water Level Predictions
The current study reveals the significance of LULC in predicting dam water levels as detected by the tested models. The assumption in the five-year time epoch used in the LULC temporal resolution is that there are insignificant changes in the natural land-covers such as water bodies, grasslands, shrublands, forests, bare soils, and land-use such as croplands. However, significant changes are expected in urban built-up, although at a slow spatial and temporal rate. Only the stochastic VAR detected the correlation and variability between the dam water levels and LULC, and predicted the dam water levels with LULC as the best predictor variable, with the highest accuracy of greater than 99%. The prediction results using MLR, RFR, and MLP-ANN showed that the LULC pattern, as interpolated over the 20-year period, may not be suitable for predicting the dam water levels for both dams as it exhibited high RMSE, MAE, and MAPE errors. For the Gaborone dam, the use of LULC resulted in a lack of convergence in prediction using the MLP-ANN. To improve the significance of LULC in dam water predictions, it is recommended to increase the temporal resolution of the LULC to annually.

Influence of Climate Factors and Climate Indices
In predicting Bokaa dam water levels using the VAR model, the combination of climate indices, rainfall, and temperature gave the best results (R 2 = 0.959, MAPE = 14.3%). This is attributed to the high correlation with the climate indices (R 2 = 0.916, MAPE = 16.2%), which resulted in good performance of all the parameters combined. Rainfall and temperature, however, did not give good results. RFR detected a higher relationship of the dam water levels using temperature (R 2 = 0.836, MAPE = 14.5%), the combination of temperature and rainfall (R 2 = 0.824, MAPE = 13.3%), the climate indices (R 2 = 0.808, MAPE = 18.3%), and the combination of climate indices, rainfall, and temperature (R 2 = 0.820, MAPE = 18.3%). While the VAR predictor variables are different for the two dams, with the exception of a combination of climate indices, rainfall and temperature, the predictor parameters for the Gaborone dam are observed to be similar to those of the Bokaa dam. It is observed that the predictions using RFR and MLP-ANN detected the variability of both dam water levels to be influenced by the same factors. For both dams using RFR and MLP-ANN, the results show that the climate factors and climate indices are the best predictors for dam water levels and are best modelled using MLP-ANN, which had the highest prediction accuracy, compared to RFR. The results further show that in the absence of reliable rainfall and temperature data, the water levels in both dams can reliably be predicted using the machine learning models based on the regional climate indices (DSLP, AI, SOI and Niño 3.4).
From the analysis of the significance of the predictor variables in Figure 15, the relatively lower contribution of rainfall in the prediction of dam water levels shows that precipitation and resulting runoff within the catchment may not be only the main sources of dam water but also marginal contributions from conjunctive water sources, such as wellfields and from other dams. As such, improvements in the prediction of the dam water levels should include the determination of the influences of the network of inter-reservoir water transfers.

Model Performances
In general, MLR was not able to detect and predict the variability of the dam water levels. On the other hand, the lower performance of VAR in detecting the influence of the seasonal climate factors and climate indices in detecting the variability of the dam water levels is attributed to the low convergence rate, as the convergence tends to be unstable, and the predictions easily fall into the local optimum trap, with an increase in the computational time, especially for the non-stationary variables [32]. On the other hand, the main advantage of the RFR machine learning, resulting in generally good results with all the variables, is in the ability to detect and discard the outlier dam water levels with ease due to the improved grouping of water level data contained in the set of terminal nodes in the decision tree. The results from MLR, VAR and RFR imply that the fluctuations in the water level in the dams are difficult to capture using the stochastic linear models [33]. The advantage of RFR and why it was able to give relatively good results is that it can handle non-linear and non-Gaussian data well and with minimal over-fitting problems as the number of trees increases [34].
MLP-ANN results support the suggestion that data-driven techniques tend to overcome the drawbacks of traditional models in terms of accuracy and the ability to model complex phenomena [35]. MLP-ANN was able to capture the influence of climate factors and climate indices with higher accuracy, though it had non-converging prediction using LULC. For the two dams, it is possible to infer that the MLP-ANN predictions adapted to the changing climate conditions. The advantages of the ANNs over other methods in predicting dam levels can be attributed to the fact that the ANN structure can detect and include the non-linear components of the system in the whole data set. Comparatively, in predicting reservoir water levels for the Angat dam in the Philippines, [25] tested the Naïve-persistence and Seasonal Mean methods as baselines against ARIMA, gradient boosting machines (GBM), and Deep Neural Networks based on LSTM, univariate (DNN-U) and multivariate models (DNN-M). The results showed that the prediction of the dam water levels was better performed using the data driven Deep Neural Network and not the traditional linear models.

VAR-ANN Hybrid Dam Water Prediction Model
The results show that neither the stochastic VAR, the decision tree based RFR, nor the MLP-ANN can independently detect the compounded impacts of LULC, climate factors, and climate indices in predicting the dam water levels. In particular, the stochastic VAR is observed to be more capable of predicting the dam water levels using LULC, which exhibits a linear trend from the five-year interval interpolations, while MLP-ANN performed better than RFR and VAR in predicting the dam water levels using the seasonal and non-linear climate factors and indices.
Since time-series hydrological data comprises different frequency components characterized by non-linear interactions, hybrid models have been proposed to improve the performance in hydrological prediction [36]. These approaches include Neural Networks based on Set Pair Analysis (SPA) and Principal Component Analysis (PCA) [37,38], Chaotic Neural Networks [39], Cluster Hybrid Neural Networks [40], And Bootstrapped Artificial Neural Networks [41,42].
From the prediction results, a hybrid dam water level prediction model comprising VAR-ANN is proposed as optimal in modeling the linear and non-linear components of the dam water levels. The VAR-ANN time-series representation of the dam water levels W L t is proposed to comprise the linear L t and non-linear N t predictor variables (Equation (11)).
N t = f e 1 , e 2 , . . . , e t−p , e t (13) In the implementation, VAR is fitted to the linear components and the outcome linearbased predictions L t at time t are derived. The residuals from the VAR, termed as ε t at time t are determined as in Equation (12). The ε t dataset after VAR fitting is considered to contain the non-linear N t time-series components of the dam water W L t levels and can be modelled using the ANN. With p input nodes, the ANN for residuals has the form in Equation (13), with f as the non-linear function estimated by the ANN and ε t is the white noise. If N t is the ANN prediction, then the hybrid prediction of at time t is defined according to Equation (14). The hybrid VAR-ANN model is implemented as depicted in Figure 17.
sidered to contain the non-linear t N time-series components of the dam water t W L levels and can be modelled using the ANN. With p input nodes, the ANN for residuals has the form in Equation (13), with f as the non-linear function estimated by the ANN and t  is the white noise. If t N  is the ANN prediction, then the hybrid prediction of at time t is defined according to Equation (14). The hybrid VAR-ANN model is implemented as depicted in Figure 17. From the best predictor variables for both dams, the average results of the hybrid VAR-ANN model for the two dams, presented in Figure 18, show an overall improvement in the prediction accuracy of the dam water levels. The results show that the hybrid model integrates the linear and non-linear variabilities in the predictor datasets to accurately predict the dam water levels. The VAR-ANN produces positive predictions using rainfall, temperature, climate indices, and LULC, with an average R 2 > 0.84 and MAPE < 10%. The results show that the average prediction RMSE, MAE, and MAPE error measures for both dams are also significantly reduced. The results imply that the hybrid model is able to capture the parametric sensitivities of both the high and low significant variables that are depicted in Figures 15 and 16. From the best predictor variables for both dams, the average results of the hybrid VAR-ANN model for the two dams, presented in Figure 18, show an overall improvement in the prediction accuracy of the dam water levels. The results show that the hybrid model integrates the linear and non-linear variabilities in the predictor datasets to accurately predict the dam water levels. The VAR-ANN produces positive predictions using rainfall, temperature, climate indices, and LULC, with an average R 2 > 0.84 and MAPE < 10%. The results show that the average prediction RMSE, MAE, and MAPE error measures for both dams are also significantly reduced. The results imply that the hybrid model is able to capture the parametric sensitivities of both the high and low significant variables that are depicted in Figures 15 and 16.

Average Model Errors and ROC Area under Curve (AUC)
The average model prediction error E (%) in Equation (15) is determined as the average for both dams using the best predictor parameters with the highest R 2 and the least RMSE, MAE, and MAPE error measures. In Figure 19, for the average predicted dam water level errors for the four models, the combination of the VAR and ANN diminishes the magnitude of the prediction error between the predicted and observed dam water levels for the two dams, producing the least errors for the predicted time-series dam water levels, and thus improved consistency in predicting the water levels.
In the first months, the E (%) for VAR-ANN is observed to be between −5% and +8% of dam water levels and diminishes to nearly 0.01% for more than 70% of the predicted dam water levels. Even though MLP-ANN performs better than RFR and VAR, its prediction errors exhibit low convergence with sinusoidal patterns in time, and this could be attributed to the influence of LULC. RFR and VAR present higher degrees of error at about 5-10%, with VAR exhibiting random spikes in error with time.
To further infer the significance of the models, the area under the receiver operating characteristic curve scores were computed for the two dams, with the results in Figure 19. The AUC scores are also based on the average true positive (sensitivity) and false positive rates (specificity) measures from the average of the best predictor variables for the dams. The results in Figure 20 show that for the Bokaa dam and Gaborone dam, VAR-ANN had the highest AUC scores, 0.89 and 0.93, performing better than MLP-ANN and RFR. The AUC scores for RFR were nearly equal, at 0.77 and 0.78, respectively, for the Bokaa dam and Gaborone dam, while VAR performance was at AUC < 0.7 for both dams. Despite the good performance from VAR-ANN, the MAPE measures for the Gaborone dam were observed to be higher than those of the Bokaa dam. The average AUC shows that the VAR-ANN has a higher ability to predict the dam water levels from all the predictor variables. Sustainability 2022, 14, x FOR PEER REVIEW 28 of 32 Figure 18. VAR-ANN model average results for dam water level predictions in Bokaa and Gaborone dams. Figure 18. VAR-ANN model average results for dam water level predictions in Bokaa and Gaborone dams. In the first months, the E (%) for VAR-ANN is observed to be between −5% and +8% of dam water levels and diminishes to nearly 0.01% for more than 70% of the predicted dam water levels. Even though MLP-ANN performs better than RFR and VAR, its prediction errors exhibit low convergence with sinusoidal patterns in time, and this could be attributed to the influence of LULC. RFR and VAR present higher degrees of error at about 5-10%, with VAR exhibiting random spikes in error with time. To further infer the significance of the models, the area under the receiver operating characteristic curve scores were computed for the two dams, with the results in Figure 19. The AUC scores are also based on the average true positive (sensitivity) and false positive rates (specificity) measures from the average of the best predictor variables for the dams. The results in Figure 20 show that for the Bokaa dam and Gaborone dam, VAR-ANN had the highest AUC scores, 0.89 and 0.93, performing better than MLP-ANN and RFR. The AUC scores for RFR were nearly equal, at 0.77 and 0.78, respectively, for the Bokaa dam and Gaborone dam, while VAR performance was at AUC < 0.7 for both dams. Despite the good performance from VAR-ANN, the MAPE measures for the Gaborone dam were observed to be higher than those of the Bokaa dam. The average AUC shows that the VAR-ANN has a higher ability to predict the dam water levels from all the predictor variables.

Conclusions
Under the influence of climate change and the intensification of land-use activities, understanding dam water capacity variations is important for planning dam water supply regimes and management. In the present study, dam water level observations in the Bokaa dam and Gaborone dam, in the semi-arid Botswana, were simulated and predicted using linear multilinear regression (MLR) and stochastic Vector AutoRegression (VAR) models, along with Random Forest Regression (RFR) and Multilayer Perceptron Neural Network (MLP-ANN) techniques. Using LULC, climate factors (rainfall and temperature) and climate indices (DSLP, Aridity Index (AI), SOI and Niño 3.4) as the dam water predictor variables, the results show that the stochastic VAR was able to detect the variation of LULC with dam water levels better than MLR, RFR and MLP-ANN, while RFR and MLP-ANN captured the relationships with the climate conditions with MLP-ANN, performing better than RFR. The stochastic VAR was not able to correlate rainfall and temperature with the dam water levels, except when integrated with the four climate indices. RFR and MLP-ANN gave the highest dam water level prediction results using rainfall, temperature, and the climate indices. MLP-ANN gave the best prediction results for the dam water level fluctuations for both dams, with the Gaborone dam predictions being more accurate than those for the Bokaa dam in terms of R 2 , but slightly lower when determined using MAPE. The higher MAPE for the Gaborone dam confirmed that the dam does not entirely

Conclusions
Under the influence of climate change and the intensification of land-use activities, understanding dam water capacity variations is important for planning dam water supply regimes and management. In the present study, dam water level observations in the Bokaa dam and Gaborone dam, in the semi-arid Botswana, were simulated and predicted using linear multilinear regression (MLR) and stochastic Vector AutoRegression (VAR) models, along with Random Forest Regression (RFR) and Multilayer Perceptron Neural Network (MLP-ANN) techniques. Using LULC, climate factors (rainfall and temperature) and climate indices (DSLP, Aridity Index (AI), SOI and Niño 3.4) as the dam water predictor variables, the results show that the stochastic VAR was able to detect the variation of LULC with dam water levels better than MLR, RFR and MLP-ANN, while RFR and MLP-ANN captured the relationships with the climate conditions with MLP-ANN, performing better than RFR. The stochastic VAR was not able to correlate rainfall and temperature with the dam water levels, except when integrated with the four climate indices. RFR and MLP-ANN gave the highest dam water level prediction results using rainfall, temperature, and the climate indices. MLP-ANN gave the best prediction results for the dam water level fluctuations for both dams, with the Gaborone dam predictions being more accurate than those for the Bokaa dam in terms of R 2 , but slightly lower when determined using MAPE. The higher MAPE for the Gaborone dam confirmed that the dam does not entirely rely on precipitation, but also on conjunctive water sources, including periodic direct supply from the Bokaa dam and wellfields. The proposed VAR-ANN hybrid model improved the prediction accuracy of the dam water levels for both dams by integrating the linear and nonlinear variabilities in the predictor datasets and the dam water levels. To improve on the current study, the temporal intervals for the LULC should be increased to annual in order to accurately capture the seasonal variabilities in the LULC; secondly, the contributions of water sources from wellfields and other dams should be incorporated into the prediction modeling. For the low convergence in the simulation and prediction of the dam water levels, using faster and hybrid tree-based machine learning algorithms is recommended for further investigations.