Spatial-Temporal Water Balance Components Estimation Using Integrated GIS-Based Wetspass-M Model in Moulouya Basin, Morocco

: The Moulouya basin in Morocco is one of many river basins around the world that are regulated with physical ﬂow control, a range of regulations and storage structures. The water budget of the basin is unbalanced; the available water resources are insufﬁcient for agricultural productivity, nature conservation and ecosystem services. This study evaluates spatial and temporal distributions of actual evapotranspiration, groundwater recharge and surface runoff for the period 2000–2020 using the WetSpass-M model in the Moulouya basin, Morocco. The WetSpass-M model’s input data are created in grid maps with the ArcGIS tool. They include monthly meteorological parameters (e.g., temperature, wind speed, rainfall,), soil map, land cover, topography, slope and groundwater depth. A good correlation has been observed between the simulated groundwater recharge and base ﬂow, with the value of R 2 = 0.98. The long-term spatial and temporal average annual precipitation of 298 mm is distributed as 45 mm (15.1%) groundwater recharge and 44 mm (14.8%) surface runoff, while 209 mm (70.1%) is lost through evapotranspiration. The simulated results showed that the average groundwater recharge of 15.1 mm (30%) falls during the summer and spring seasons, while the remaining 29.5 mm (70%) occurs during the winter and autumn seasons. Annually, 2430 million m 3 of water recharges to the groundwater system from the rainfall for the entire basin. The study’s ﬁndings would help local stakeholders and policymakers in developing sustainable and effective management of available surface water and groundwater resources in the Moulouya basin.


Introduction
Groundwater is an essential supply of fresh water for drinking, industry, agriculture, and the preservation of ecosystem services all over the world [1,2]. Human intervention alters the hydrological process at increasing rates which, in return, reduces the availability of water. Furthermore, climate change, rising global population and shortages of groundwater recharge all contribute to lowering groundwater levels. A better understanding of the special and temporal distributions of water balance components, particularly groundwater recharge, is important for successful water resource management and modeling subsurface fluid and contaminant transport [3], especially as these resources become the primary source of drinking water [4].
As a semi-arid country, Morocco is characterized by erratic precipitation. Hence, groundwater and surface water resources are essential for socio-economic growth [5]. As a result, aquifers are extensively exploited to meet growing agricultural, domestic, industrial and tourist demands [6][7][8]. Unsustainable land-use practices such as intensive farming, over-exploitation of forest resources and overgrazing, particularly by the poor communities, are significantly altering the hydrologic features of the Moulouya basin.

Hydrological Modelling (WetSpass-M)
Estimating groundwater recharge is critical for sustainable and effective management of surface water and groundwater resources. Estimation of groundwater recharge is a challenging and complex process because groundwater recharge relies on several parameters such as topography, soil texture, land use/land cover (LULC), groundwater depth, meteorological parameters and other hydrologic characteristics [41]. The WetSpass model was developed as a physically-based approach for estimating the long-term average actual evapotranspiration, groundwater recharge and surface runoff [22,41]. In this study, a WetSpass-M model [42] is applied to evaluate the spatial distribution of water balance components on monthly, seasonal and annual scales. The model takes into account the spatial distribution of elevation, slope soil texture, LULC and meteorological parameters for each raster cell. The WetSpass model treats a region or basin as a consistent pattern of grid cells [41]. The water balance parameters of each area is obtained by the vegetated, open-water, bare soil, while impervious fractions per grid pixel are determined by the following equation [22]: pattern of grid cells [41]. The water balance parameters of each area is obtained by the vegetated, open-water, bare soil, while impervious fractions per grid pixel are determined by the following equation [22]:  Figure 2. The detailed description of WetSpass-M is presented in appendix A.

Input D,iata
The required input data for the WetSpass-M model can be categorized into two groups: parameters tables and GIS grid maps [22]. GIS grid maps include meteorological data (rainfall, potential evapotranspiration (PET), average temperature and wind speed), groundwater depth, slope, topography, soil and LULC type. The LULC and soil type are linked to the model using attribute tables of the soil and LULC raster maps. Moreover, the attribute tables make it simple to define new soil or LULC types, as well as to adjust parameter values, allowing for the study of future water and land management scenarios [41]. The input data were prepared as a raster map in ESRI ASCII grid format based on the digital elevation model (DEM) with a cell size of 100 × 100 m totaling (3886,3306) raster cells. Table 1 presents the input parameters for the WetSpass-M model. Wind speed MHBA and own processing 100 × 100 m

Input Data
The required input data for the WetSpass-M model can be categorized into two groups: parameters tables and GIS grid maps [22]. GIS grid maps include meteorological data (rainfall, potential evapotranspiration (PET), average temperature and wind speed), groundwater depth, slope, topography, soil and LULC type. The LULC and soil type are linked to the model using attribute tables of the soil and LULC raster maps. Moreover, the attribute tables make it simple to define new soil or LULC types, as well as to adjust parameter values, allowing for the study of future water and land management scenarios [41]. The input data were prepared as a raster map in ESRI ASCII grid format based on the digital elevation model (DEM) with a cell size of 100 × 100 m totaling (3886,3306) raster cells. Table 1 presents the input parameters for the WetSpass-M model. The meteorological data for the period 2000-2020 were obtained from the Moulouya Hydraulic Basin Agency. The temperature varies greatly depending on the altitude, the orientation of the valleys or basins, or the geographical location within the mountain range. The climate is classed as semi-arid: moderate temperate with winter precipitation [43]. Moulouya basin is one of the driest regions in Morocco with an inconsistent annual rainfall of 330 mm [39,44]. The average monthly temperature of the basin is 16 • C, ranging from 7 • C to 20 • C (Figure 3a). January is the coldest month, while July represents the warmest month. Low temperatures are observed on the High Atlas ridges in the South-West and Middle Atlas Mountains in the West, while higher temperatures are concentrated along the Mediterranean coast and passage Missour-Outat El Hadj. Effective precipitation is a major factor in the direct recharge of groundwater through infiltration and vertical percolation through the saturated zone [45]. The average monthly precipitation of the Moulouya basin ranges between 2.1 mm and 36.7 mm with an average value of 25 mm. Average annual rainfall ranges from 154 mm yr −1 to 770 mm yr −1 , with a mean rate of 298 mm yr −1 (Figure 3b). About 60% of the rainfall occurs in the autumn and winter seasons, while the remaining 40% takes place in the spring and summer seasons. The highest rainfall is observed in the mountain Middle Atlas (Bounacer) and high atlas (Ayachi), while the lowest values are recorded in Missour Outat El Hadj Corridor and the High Plateaus.
nel2A Land Cover (ESA-S2-LC20) 7 Soil texture FAO database 100 × 100 m 8 Groundwater depth MHBA and own processing 100 × 100 m 9 lookup The meteorological data for the period 2000-2020 were obtained from the Moulouya Hydraulic Basin Agency. The temperature varies greatly depending on the altitude, the orientation of the valleys or basins, or the geographical location within the mountain range. The climate is classed as semi-arid: moderate temperate with winter precipitation [43]. Moulouya basin is one of the driest regions in Morocco with an inconsistent annual rainfall of 330 mm [39,44]. The average monthly temperature of the basin is 16 °C, ranging from 7 °C to 20 °C (Figure 3a). January is the coldest month, while July represents the warmest month. Low temperatures are observed on the High Atlas ridges in the South-West and Middle Atlas Mountains in the West, while higher temperatures are concentrated along the Mediterranean coast and passage Missour-Outat El Hadj. Effective precipitation is a major factor in the direct recharge of groundwater through infiltration and vertical percolation through the saturated zone [45]. The average monthly precipitation of the Moulouya basin ranges between 2.1 mm and 36.7 mm with an average value of 25 mm. Average annual rainfall ranges from 154 mm yr −1 to 770 mm yr −1 , with a mean rate of 298 mm yr −1 (Figure 3b). About 60% of the rainfall occurs in the autumn and winter seasons, while the remaining 40% takes place in the spring and summer seasons. The highest rainfall is observed in the mountain Middle Atlas (Bounacer) and high atlas (Ayachi), while the lowest values are recorded in Missour Outat El Hadj Corridor and the High Plateaus. Evapotranspiration is one of the most crucial factors in the water cycle but also one of the most spatially distributed variables [46]. The evapotranspiration values available on a monthly scale have also been downloaded from the Global Land Data Assimilation System (GLDAS) and Famine Early Warning Systems Network (FEWS NET) Land Data Assimilation System (FLDAS) data set. The average annual evapotranspiration ranges from 150 mm to 500 mm, with a mean value of 207.22 mm. The average annual wind speed Evapotranspiration is one of the most crucial factors in the water cycle but also one of the most spatially distributed variables [46]. The evapotranspiration values available on a monthly scale have also been downloaded from the Global Land Data Assimilation System (GLDAS) and Famine Early Warning Systems Network (FEWS NET) Land Data Assimilation System (FLDAS) data set. The average annual evapotranspiration ranges from 150 mm to 500 mm, with a mean value of 207.22 mm. The average annual wind speed of the Moulouya basin varies from 3.6 m/s as a minimum and 6.4 m/s maximum values with an average speed of 5 m/s (Figure 3c). Monthly groundwater data for 189 wells are obtained from Moulouya Hydraulic Basin Agency (MHBA), for the period (2000 to 2020). Kriging interpolation was used to develop the spatial distribution of average monthly groundwater depth ( Figure 3d). The groundwater depth ranges between 10.6 m to 80.7m with an average depth of 30.5 m.

Topography and Slope
Most of the literature report that geomorphology is the most important factor for groundwater [47][48][49]. The DEM of the Moulouya basin is obtained from SRTM ( Figure 4a). The highest point of the investigated area is 3730 m in the northern part of the Gurage Mountains while the lowest is 2 m in the southern part of the basin. The average elevation of the investigated basin is 1140 m. The slope gradient directly controls the surface water infiltration. A low slope gradient restricts water flow and hence increases infiltration rate [50], whereas a high slope gradient has limited groundwater recharge owing to high surface runoff [51]. Slope analysis tool in the ArcGIS environment was applied to create the slope map from the DEM. The slope ranges from 0 to 75 degrees with an average degree of 5.9 ( Figure 4b).

Land Use
LULC type has a substantial effect on groundwater recharge or infiltration [52,53]. The LULC is also useful for estimating the values of the vegetative parameters such as evaporative zone depth and leaf area index (LAI). The LAI parameter is used to control

Land Use
LULC type has a substantial effect on groundwater recharge or infiltration [52,53]. The LULC is also useful for estimating the values of the vegetative parameters such as evaporative zone depth and leaf area index (LAI). The LAI parameter is used to control both surface evaporation and transpiration [54]. LULC is one of the most essential controlling factors in basin hydrology [55]. LULC data for the Moulouya basin are obtained from ESA-S2-LC20 with a resolution of 20 m (http://2016africalandcover20m.esrin.esa.int/download. php, accessed on 10 April 2021) (Figure 4c) [56]. This classification, which identified 10 generic classes ("shrubs", "trees", "vegetation aquatic" "grassland", "cropland", "sparse vegetation/mosses and lichen", "built-up", "bare areas", "snow", and "open water"), was the result of combining the machine learning (ML) and random forest algorithms. Eight land-use classes are identified in the investigated area, which are reclassified into seven land-use classes for the WetSpass database. Seven classes of land cover have been identified the investigated basin as shown in (Figure 4c). It is dominated by a shrub area of (35.16%), a bare surface of (34.11%), a heather surface area of (24.86%), agriculture areas (4.98%) and a total area of mixed forest, build up and water bodies of (0.9%). In recent years, the upper Moulouya Basin has witnessed a shift toward major irrigated areas where the cultivation of fruit trees (apples) is concentrated. In addition to the remarkable development of the olive, tree lands have been observed in the Central Moulouya Basin, which relies primarily on groundwater for irrigation.

Soil Data
Soil texture is the key point for understanding all required information for the hydrological investigation of any region [57]. Soil represents the main physical characteristics that control runoff and recharge. Soil infiltration capacity relies on soil permeability, which calculates its storage capacity and influences the hostility of the flowing water into deep layers. Soil texture has a significant impact on the soil infiltration capacity. Sandy soil has the highest infiltration rates while heavy clay and loamy soil shows a lower infiltration rate and higher surface runoff [58]. The soil map of the Moulouya basin is obtained from the Harmonized World Soil Database (HWSD) [59] (Figure 4d). The dominant soil texture of the basin is clay loam which covers 92.79% (48,025.10 Km 2 ) of the total study area, while loam and clay represent 7.10% (3669.6 Km 2 ) and 0.12% (58.83 Km 2 ) of the basin, respectively.

Validation of WetSpass-M Model
The Moulouya basin is subdivided into 23 sub-catchments. Within these regional catchments ( Figure 5), the daily discharge of 17 river gauging stations for the period 2000-2020 is obtained from MHBA to perform the hydrograph analysis. Based on the DEM and river network map of the Moulouya basin, we utilize the HydrologyTool in GIS to delineate the areas associated with each gauging station. The automated Web-Based Hydrograph Analysis Application (WHAT) is applied to derive a base flow from streamflow data. WHAT has three separating filters: the Eckhardt recursive digital filter method (RDF) [60], the one-parameter digital filter method (OPM) [61][62][63][64] and the localminimum method (LMM) [65].
The Eckhardt recursive digital filter method (RDF) [60] is applied in this study: where b t represents base flow at time step t (m 3 /s); b t−1 represents the filtered base flow at time step t−1 (m 3 /s); BFI max presents the maximum long term ratio of base flow/total streamflow; Q t is the total streamflow at time step t (m 3 /s) and α is the filter parameter. Eckhardt [60] suggested BFImax values of 0.50 for ephemeral streams including porous aquifers, 0.25 for perennial streams containing hard rock aquifers and 0.80 for perennial streams containing porous aquifers. These results were acquired by using and validating this filtering method on watersheds in Maryland, Germany, Illinois and Pennsylvania [65].
The proposed values of 0.80 for BFImax and 0.98 for the filter parameter, which correlate to the hydrogeological characteristics of the region, were used in this case.

Validation of WetSpass-M Model
The calculated base flow rates with RDF method for daily streamflow records at 17 gauging stations between 2000 and 2020 are used to validate the WetSpass-M model. The simulated groundwater recharge for the WetSpass-M model at the corresponding gauging stations is extracted from the spatially distributed results in GIS. Figure 6 demonstrates that the simulated average monthly groundwater recharge by WetSpass-M matches the base flow with R 2 = 0.98, a mean error of 7 mm month −1 and an absolute mean error of 18 mm month −1 . According to the findings, the agreement between simulated and measured recharge lies within an acceptable range.

Water Balance Components
The WetSpass-M model's major outputs are ASCI maps of monthly real evapotranspiration, surface runoff and groundwater recharge, from 2000 to 2020 (252-time steps). Every pixel on these maps shows the water balance component's magnitude (in mm). This is the first investigation of water balance components in the Moulouya basin. An assessment of annual water balance components is required to assess the water budget of the Moulouya basin, also for seasonal and monthly scales to calculate the water needs for agriculture. The WetSpass findings for water balance components will be utilized as inputs and boundary conditions to integrate groundwater modeling of the Moulouya basin. The simulated annual, seasonal and monthly water balance components are shown in Table 2.  The total actual evapotranspiration per pixel is determined by a WetSpass-M model as the sum of evaporations from bare soil, impermeable surface area, open water, interception and transpiration of the vegetated area [15,42]. The WetSpass-M estimates the long-term monthly actual evapotranspiration of the Moulouya basin to be 6.5 mm and 78.6 mm, as the lowest and highest values, respectively, with a mean value of 17 mm month −1 and a standard deviation of 4 mm month −1 ( Table 2). Annual actual evapotranspiration is determined by summing the monthly values over the year. The minimum, maximum and mean annual actual evapotranspiration for the studied period are 78 mm, 943 mm and 209 mm, respectively. Mean annual actual evapotranspiration attributes to 70.1 % of the average annual precipitation (Figure 7a). Mean long-term actual evapotranspiration values for the wet (autumn and winter) seasons and dry (spring and summer) seasons are 130.5 mm, 87.9 mm, respectively. About 60% of the annual actual evapotranspiration occurs in the autumn and winter seasons while the other 40% takes place in the spring and summer seasons. This difference is due to rainfall variation within the two seasons. Actual evapotranspiration ranges from 1.9 mm to 27.4 mm in July and April, respectively, with a mean monthly value of 17.4 mm month −1 (Figure 8). The percentage of actual evapotranspiration to rainfall in the basin is highest for July (90.4%) and June (80.2%). As depicted in Figure 7a, high seasonal and annual actual evapotranspiration are recorded in the southwest of the Moulouya because of higher precipitation, while the northeast part that receives less rainfall shows lower evapotranspiration.  The spatial variability of average annual surface runoff is given in Figure 7b. Such a monthly average surface runoff varies from 0 mm yr −1 to 53 mm yr-1, with a mean rate of 4 mm yr-1 and a standard deviation of 13 mm month −1 ( Table 2). The average annual surface runoff is calculated depending on simulated monthly data. Annual surface runoff presents a large spatial variation with rates between 5 mm and 635 mm. The mean and standard deviation of this variability is 44 mm yr −1 and 40 mm y −1 , respectively ( Table 2). Average annual surface runoff accounts for 14.9% of the total average annual rainfall (Figure 7b). About 80% of the simulated surface runoff takes place in wet seasons (winter and autumn) while the other 20% occurs in dry seasons (spring and summer). Surface runoff varies from 0.1 mm to 7.5 mm in July and December, respectively. Surface runoff is highest in December, which is the month with the highest rainfall while the lowest is recorded in May, June and July, which coincide with the month with the lowest rainfall ( Figure 8). The highest mean annual and seasonal surface runoff values are recorded in the central west of the Moulouya basin which is attributed to their gentle slope and the presence of loam, clay loam and clay soils with low permeability (Figure 7b). The spatial variability of average annual surface runoff is given in Figure 7b. Such a monthly average surface runoff varies from 0 mm yr −1 to 53 mm yr −1 , with a mean rate of 4 mm yr −1 and a standard deviation of 13 mm month −1 ( Table 2). The average annual surface runoff is calculated depending on simulated monthly data. Annual surface runoff presents a large spatial variation with rates between 5 mm and 635 mm. The mean and standard deviation of this variability is 44 mm yr −1 and 40 mm y −1 , respectively ( Table 2). Average annual surface runoff accounts for 14.9% of the total average annual rainfall (Figure 7b). About 80% of the simulated surface runoff takes place in wet seasons (winter and autumn) while the other 20% occurs in dry seasons (spring and summer). Surface runoff varies from 0.1 mm to 7.5 mm in July and December, respectively. Surface runoff is highest in December, which is the month with the highest rainfall while the lowest is recorded in May, June and July, which coincide with the month with the lowest rainfall ( Figure 8). The highest mean annual and seasonal surface runoff values are recorded in the central west of the Moulouya basin which is attributed to their gentle slope and the presence of loam, clay loam and clay soils with low permeability (Figure 7b). Groundwater recharge is crucial for evaluating groundwater resources; nevertheless, it is complex to estimate recharge [66,67]. The spatial distribution of groundwater recharge depends on meteorological conditions, topography, slope, groundwater depth, land cover and soil type. The WetSpass-M model simulates groundwater recharge for the Moulouya as a residual parameter of the water balance parameters by subtracting the evapotranspiration and runoff from the monthly precipitation. Seasonal groundwater recharge of the Moulouya basin varies spatially with the topography and characteristics of the basin (Figure 9a-d). The simulated monthly groundwater recharge of the Moulouya basin ranges from 0 mm month −1 to 29 mm month −1 . The mean and standard deviations are 4 mm and 2 mm, respectively. The average annual groundwater recharge is calculated based on monthly simulated data. The annual average of groundwater recharge varies from 0 mm yr −1 to 350 mm yr −1 , with a mean value of 45 mm yr −1 and a standard deviation of 22 mm yr −1 ( Table 2). Average groundwater recharge represents 15% of average annual rainfall (Figure 7c), of which an average of 15.1 mm (30%) falls during spring and summer seasons, while the other 29.5 mm (70%) occurs in the winter and autumn seasons ( Figure 9 and Table 2). The estimated monthly groundwater recharge in the Moulouya basin is presented in Figure 8. The average monthly groundwater recharge ranges from 0.1 mm in July to 16.6 mm in January. The highest average rainfall amount occurs in January, of which 18.5% infiltrates into the groundwater system ( Table 2). As shown in Figures 7c and 10, the mountainous areas in the central-western part of the Moulouya basin (Medium Atlas, High Atlas, Benisnasen), that receive a high amount of rainfall has higher seasonal and annual and groundwater recharges. The highest groundwater recharge is detected in forest and agricultural areas in the central and southern parts of the Moulouya basin. On the contrary, the northern part accounts for lower groundwater recharge which is attributed to the presence of theater and bare areas with less-permeable clay soils (Figure 9a-d). Figures 10 and 11 show the water balance parameters (actual evapotranspiration, groundwater recharge and surface runoff) as a function of various LULC types in soil textures, respectively. Actual evapotranspiration, groundwater recharge and surface runoff depend on LULC. Forest areas are characterized by high groundwater recharge with an average of 208 mm yr −1 , while the heather areas have the lowest surface runoff with an average of 26.4 mm yr −1 (Figure 10). Open water is assigned by a zero-groundwater recharge in the WetSpass model, assuming that the groundwater recharge obtained from precipitation on the open water portion is insignificant in comparison to the probable recharge from the surface water [41]. Shrub and heather areas show high average actual evapotranspiration by 229 and 225 mm yr −1 , respectively, while they have average groundwater recharges of 49.67 and 48.25 mm yr −1 , respectively ( Figure 10). Built-up areas have a low groundwater recharge and actual evapotranspiration with an average of 45 mm yr −1 , respectively, as these areas are characterized by a partially or fully impervious surface. Water balance components are also strongly dependent on soil textures. Spatial variability of soil textures has a considerable impact on local/regional hydraulic characteristics [68]. Textural characteristics influence hydraulic conductivity, specific yield, porosity and capillarity. The Moulouya basin is characterized by heavy soils (loam, clay loam, clay). Heavy soils (clay and clay loam) have the highest surface runoff and evapotranspiration. The groundwater recharge value of clay is approximately two-thirds that of loamy soil textures. The WetSpass-M model simulates the average annual actual evapotranspiration and surface runoff of clay soils to be 219 and 55 mm yr −1 , respectively, while the simulated average annual groundwater recharge of clay soils is 40 mm yr −1 (Figure 11).

Water Balance Components under Different Land Use/Land Cover (LULC) Uses and Soil Textures
This work shows that water balance components have been affected by both LULC classes soil texture Tables 3, A1 and A2. Table A1 shows the average annual evapotranspiration (mm) across various combinations of LULC classes and soil types. The highest mean annual evapotranspiration is found in regions dominated by water bodies, heather and shrub, whereas areas with built-up, forest and farming activities have the lowest evapotranspiration rate under all soil textures (Table A1). This could be due to the high transpiration demand of vegetation cover and water availability of soil type. The current findings also reveal that land-use variability affects the spatial variability of actual evapotranspiration more than soil type (Table A1). The smaller variability in average yearly evapotranspiration among soil types shows that the evapotranspiration rate in the studied basin is less reliant on soil texture. Singh et al. [69] report similar results of the impacts of soil texture and LULC classes on evapotranspiration variability.  Higher surface runoff is observed from the built-up and bare areas (Table A2). This is because heavy soil has a lower infiltration capacity and rural settlement land units are more compact, resulting in a decreased rate of water recharge and percolation [70,71]. On the other hand, areas dominated by loamy soil with heather and shrub coverage show a lower surface runoff (Table A2). Similar results on the impacts of vegetation and soil texture on spatial variability of runoff are reported by [72,73]. Consequently, the spatial and temporal distribution of surface runoff can help in understanding the primary parameters that drive runoff variability in Moulouya basin conditions. Moreover, the present findings indicate that LULC variability affects the spatial distribution of surface runoff more than the impact of soil type in the basin. Similar results of the impacts of LULC classes and soil textures on surface runoff variability are reported by [74][75][76]. Table 3 shows the average annual groundwater recharge (mm) across different combinations of soil types and LULC classes to evaluate the spatial variations of the groundwater recharge under different soil textures and LULC classes in the Moulouya basin. Forest and agriculture areas with loamy soils have a higher groundwater recharge contributing to the high permeability of loam textured soil compared to clay soil. The highest average groundwater recharge of 220.6 mm is accounted for by forest areas with loamy soil type while bare areas with clay soils show the lowest average groundwater recharge of 21.3 mm followed by heather areas ( Table 3). The results indicate that groundwater recharge is more affected by LULC than soil type since the standard deviation of the groundwater recharge for the various LULC classes is higher than the standard deviation of groundwater recharge for the various soil textures (Table 3). Similar results of groundwater recharge variability with LULC classes and soil textures to the present finding are reported for the Flanders (Belgium) by [77], for the upper San Pedro watershed, Mexico and USA, by Nie et al. [78].

Conclusions
The Moulouya basin is significantly affected by human activities, and it suffers from water stress and intensification of drought hazards. Developing a groundwater model for the basin requires an accurate evaluation of evapotranspiration and groundwater recharge as boundary conditions. A GIS-based water balance model, WetSpass-M, was applied to estimate annual seasonal and monthly groundwater recharge, actual evapotranspiration and surface runoff in the Moulouya basin for the period from 2000 to 2020. The total water balance components of the bare soil, vegetated, open-water and impervious fraction per grid pixel were determined. The WetSpass-M model's main input variables included climate data (air temperature, precipitation, potential evapotranspiration and wind speed), LAI, soil types, groundwater depth, DEM, slope and LULC of the investigated area. Such input data were created as raster maps using ArcGIS environment. The spatial variability of recharge depends on topography, slope, climate conditions, groundwater depth, LULC and soil texture. LULC and soil textures were found to be dominated by agricultural areas and clay-loamy soils in the Moulouya basin. Water balance components under various LULC and soil textures were evaluated. The base flow of the Eckhardt recursive digital filter method at 17 gauging stations was implemented to validate the efficiency of the WetSpass-M.
The WetSpass-M model evaluated the annual surface runoff of the basin, for the period from 2000 to 2020, to be 5 mm and 635 mm as minimum and maximum values respectively, which attributes for 14.9% of the mean annual precipitation. Around 15% (199 mm yr −1 ) of the average annual precipitation is accounted for by groundwater recharge with minimum and maximum recharge of 78 mm yr −1 and 943 mm yr −1 , respectively. Annually simulated actual evapotranspiration varies from 175 mm yr −1 to 412 mm yr −1 with an average of 209 mm yr −1 . This represents 70.1% of the annual average precipitation; 80% of total actual evapotranspiration takes place in the wet season, while the other 20% occurs during dry seasons. Simulation outputs indicate that the WetSpass-M model was used correctly to estimate water budget components in the Moulouya basin. This study can be utilized to develop an integrated groundwater modeling and to evaluate possible areas for controlled artificial recharge by runoff harvesting to increase groundwater storage. Artificial recharge may be used to retain stormwater runoff during rainy seasons to increase groundwater availability during dry months. Additionally, monthly groundwater withdrawals should be managed with respect to recharge patterns, both spatial and temporal.
where I m is the interception (mm/month), P m is monthly rainfall [mm/month] and I R is interception ratio. Monthly interception ratio, I R is determined according to [79]: where d P is the number of rainy days per month and I D the daily interception threshold, which depends on land use [80].
where LAI is the leaf area index and a is an interception parameter.

Appendix A.2. Surface Runoff
The monthly surface runoff SRm was calculated using a rational method using two coefficients [41]: where C sr actual runoff coefficient [-], C h is a coefficient [-] representing soil moisture condition [81].
θ s represents soil moisture content [m 3 /m 3 ], θ sat is the soil porosity [m 3 /m 3 ] and b is an exponent [-] representing the effect of precipitation intensity.

Appendix A.3. Evapotranspiration
Evapotranspiration per raster cell is estimated by adding up actual evapotranspiration from the bare soil (ET b ) vegetated (ET v ), open water (Eo) and impervious surface (ET i ). WetSpass-M uses potential evaporation and vegetation coefficients for calculating the actual evapotranspiration. Vegetation coefficient is required to estimate the reference transpiration from the potential evapotranspiration (ET P ): In vegetated areas, the actual transpiration is modified when the groundwater table is below the root zone where a i is a calibrated parameter (adopted from WetSpass) related to the sand content of the soil texture and w is the available water for transpiration, which is given by: where R d is the rooting depth, θ f c − θ pwp is the available water content in the plant per time step, represented as the difference in water content at field capacity and permanent wilting point. Summation of actual transpiration and interception produces the actual evapotranspiration for the vegetated area (ET v ). The total monthly evapotranspiration (ET m ) [mm/month]) per cell is: with, respectively, area fraction and evapotranspiration for bare soil a b , ET b ; vegetated area a v , ET V , open water a o , ET O ; and impervious surface a i , ET i [22,41].

Appendix A.4. Recharge
The monthly recharge (R m (mm/month)) WetSpass-M is determined as the residual term of the monthly water balance: Appendix B