A Remote Sensing-Based Assessment of Water Resources in the Arabian Peninsula

: A better understanding of the spatiotemporal distribution of water resources is crucial for the sustainable development of hyper-arid regions. Here, we focus on the Arabian Peninsula (AP) and use remotely sensed data to (i) analyze the local climatology of total water storage (TWS), precipitation, and soil moisture; (ii) characterize their temporal variability and spatial distribution; and (iii) infer recent trends and change points within their time series. Remote sensing data for TWS, precipitation, and soil moisture are obtained from the Gravity Recovery and Climate Experiment (GRACE), the Tropical Rainfall Measuring Mission (TRMM), and the Advanced Microwave Scanning Radiometer for Earth Observing System (AMSR-E), respectively. The study relies on trend analysis, the modiﬁed Mann–Kendall test, and change point detection statistics. We ﬁrst derive 10-year (2002–2011) seasonal averages from each of the datasets and intercompare their spatial organization. In the absence of large-scale in situ data, we then compare trends from GRACE TWS retrievals to in situ groundwater observations locally over the subdomain of the United Arab Emirates (UAE). TWS anomalies vary between − 6.2 to 3.2 cm/month and − 6.8 to − 0.3 cm/month during the winter and summer periods, respectively. Trend analysis shows decreasing precipitation trends ( − 2.3 × 10 − 4 mm/day) spatially aligned with decreasing soil moisture trends ( − 1.5 × 10 − 4 g/cm 3 /month) over the southern part of the AP, whereas the highest decreasing TWS trends ( − 8.6 × 10 − 2 cm/month) are recorded over areas of excessive groundwater extraction in the northern AP. Interestingly, change point detection reveals increasing precipitation trends pre- and post-change point breaks over the entire AP region. Significant spatial dependencies are observed between TRMM and GRACE change points, particularly over Yemen during 2010, revealing the dominant impact of climatic changes on TWS depletion. observed after change point division of the TRMM time series. Change point detection revealed signiﬁcant spatial dependencies between the TRMM and GRACE products, especially around Yemen during 2010, indicating hotspots of climate change impacts on TWS variations that should be accounted for alongside anthropogenic factors.


Introduction
The Arabian Peninsula (AP) has been experiencing a significant increase in water demand as a result of population growth, urbanization, and an overall increase in living standards. The stress on the region's already scarce water resources is further exacerbated by the increasing reliance on groundwater, especially in deep inland locations where water supply from seawater desalination is not economically feasible. Even in coastal areas, the frequent deterioration of seawater quality due to harmful algal blooms or oil spills [1,2] hinders the operation of desalinations plants. In such context, the monitoring of renewable water resources is essential, yet also challenging, due to the sporadic nature of rainfall events that are extremely unevenly distributed spatially and temporally over the AP [3][4][5][6][7]. Studying small watersheds within the AP cannot reflect the macroscale processes and capture the flow of the dominant groundwater component, which is important in arid regions given the absence of sustained surface water flow [8]. In turn, examining the variability over a regional scale is essential to understand the interaction between key hydrological variables and their dynamics in the region.
The spatial distribution of the selected key hydrological variables is also investigated and the temporal consistency of the detected changes in their variability is analyzed. In the absence of large-scale and dense observation networks over the entire study area, the trend of GRACE estimates is verified locally using observations from wells in the United Arab Emirates (UAE). The inferred agreement between GRACE and local well observations in the UAE is assumed to be reflective of the performance of the retrieval of TWS anomalies over the entire study area. Moreover, here, we test the feasibility of inferring relationships between key hydrological variables and macroscale processes by exclusively using remote sensing data. This makes the proposed analysis expandable to other ungauged watersheds where in situ observations are scarce or not available. Figure 1a shows the topography of the AP varying from the lowest points along the coast of the Arabian Gulf to over 3500 m at the highest southern peaks in Yemen. The Hajar Mountain range on the border between Oman and the UAE is also depicted with peaks reaching 3000 m. Recharge mechanisms are highly dependent on both the existing topography and soil textures. Figure 1b shows widespread shallow and poorly developed soils, enriched with lime, gypsum, or salts, which are strongly coupled to the arid desert climate. Sand dunes and sheets of transported materials cover the majority of the region. Alternatively, shallow and stony soils cover the southern highlands in Yemen and the Hajar Mountain range, which may favor lateral flow rather than recharge of the deep aquifer. Ephemeral rivers (wadis) are fed by runoff from these highlands, which frequently trigger flashfloods downstream [25][26][27]. The absence of a dense network of gauged streams hinders the verification of hydrological processes, which makes the reliance on remotely sensed products an attractive alternative for regional studies over the AP.

Study Area
The dry (June, July, August, and September-JJAS) and wet (December, January, February, and March-DJFM) seasons for much of the AP are consistent with the region and typical hydrological regimes. However, accurate details on the seasonal rainfall distribution for the AP using long-term observations are generally not available in the literature [5]. Large rainfall amounts are observed over the northern AP during the wet season associated with pressure troughs from the eastern Mediterranean, while the southern portion records peak rainfall amounts during the dry season with the advancement of moist air masses from the Indian Monsoon Trough [30][31][32]. Consequently, two transitional periods (April-May and October-November) separate these dry and wet regimes, which dictate the region's climatology. Remote Sens. 2021, 13, x FOR PEER REVIEW 4 of 18 Figure 1. Terrain elevations (a) and soil types (b) from the GTOPO30 Digital Elevation Model [28] and De Pauw [29], respectively. The spatial distribution of the observational well network in the United Arab Emirates (UAE) is shown in subplot (c).  [28] and De Pauw [29], respectively. The spatial distribution of the observational well network in the United Arab Emirates (UAE) is shown in subplot (c).

Remote Sensing Products
The TRMM 3B42 Research Derived Daily Product is used for precipitation data over the time period extending from January 1999 to July 2015 [33]. The algorithm applied for this product is the Version 7 TRMM Multi-Satellite Precipitation Analysis (TMPA), which is a major improvement over the previous version 6, especially in terms of daily resolution [34]. This product has been shown to capture light rain events and outperform other global precipitation products when intercompared over the AP and surrounding areas [3,9,11]. For soil moisture, level 3 daily AMSR-E estimates are obtained over the period of June 2002 to September 2011 [35,36]. AMSR-E estimates of soil moisture have been found to perform well over desert environments in the absence of vegetation [37,38] and were applied in conjunction with the TRMM product for flood monitoring over Saudi Arabia [39]. Both TRMM and AMSR-E products have a high spatial resolution of 0.25 • × 0.25 • (latitude × longitude) at the daily scale. Three GRACE solutions using the standard spherical harmonic approach (SSH) [40] are released from different processing centers, namely, GeoForschungsZentrum Potsdam (GFZ), Center for Space Research at University of Austin Texas (CSR), and Jet Propulsion Laboratory (JPL). While slight differences have been reported between the three products, their trends remain consistent [41]. Here, we use the JPL GRACE product with scaling coefficients obtained from the National Center for Atmospheric Research's numerical land-hydrology community model (NCAR CLM4). The GRACE data time period extends from April 2002 to October 2017, with a resolution of 1 • × 1 • at the monthly scale [40,42]. For the sake of the intercomparison, the products are analyzed for the overlapping period of June 2002 to September 2011. Table 1 includes a summary of the used datasets.

In-Situ Well Observations
In situ well observations collected by the UAE Ministry of Climate Change and Environment (MoCCAE), formerly known as the Ministry of Environment and Water (MOEW), are used to locally test the inferred trends from GRACE. Given the lack of such observations across the entire AP, an agreement between well observations and GRACE anomalies over the test site in the UAE should indicate the validity of expanding the GRACE-based analysis to the rest of the study area. The records include monthly water table depths of 39 wells logged from January 1989 to September 2014. Wells in the UAE are located around a network of small dams and flood wave breakers built on major wadis to maximize groundwater recharge. The elevation of each well, with respect to sea level, is also recorded. After subtracting the well elevations, monthly water table elevations for each well are determined, starting from April 2002 to match the GRACE temporal coverage. Next, a monthly spatial average of the water table variations across the wells is computed, from which average water storage anomalies are derived using the same storage reference period as GRACE, namely, January 2004 to December 2009.

Data Processing
A land mask of the AP is applied to the processed datasets. A spatial analysis of the derived climatology from the GRACE, TRMM, and ASMR-E products from 2002 to 2011 for winter (DJFM) and summer (JJAS) is first carried out by examining maps of their spatial patterns to detect possible similarities. Then, the TRMM and AMSR-E products are aggregated to match the temporal resolution of the GRACE data resolution by calculating the cumulative monthly sum on a pixel-per-pixel basis. Trend analysis is carried out using the native spatial and temporal resolutions with the exception of the AMSR-E product, for which the aggregated monthly dataset is used to avoid potential errors associated with orbital gap filling. The final trend coefficients are mapped after interpolation by kriging to express the spatial variation in terms of the variogram and minimize prediction errors [43].

Trend Analysis
The time series of the three datasets are used to fit long-term linear trends to the evolution of water storage, precipitation and soil moisture over the AP. The derived trends represent the overall 10-year variation through both the magnitude and the sign of the resulting trend coefficients (slopes) at each pixel across the study domain. Next, the modified Mann-Kendall [44] test is carried out on a pixel-per-pixel basis for the three datasets to test for their trend significance (5% significance level) throughout the study period. This test accounts for autocorrelations present in the data by modifying the variance of the non-parametric Mann Kendall Test [45,46]. For a given sample x 1 , x 2 , . . . , x n , of size n, the standardized test statistic Zs is obtained as: and with m being the number of tied values in t i. The null hypothesis of trend absence is rejected at a significance level α (5% used here) if the p-value is less than α.

Change Point Analysis
The cumulative sum (CUSUM) method is a simple nonparametric statistical technique that can be used for detecting abrupt change points in time series [47]. This method has been shown to be computationally inexpensive while remaining statistically optimal [48]. For a given time series, x 1 , x 2, . . . , x n , the cumulative sum of deviations at any time k is given by: Change point detection using the CUSUM method is carried out on a pixel-per-pixel basis across the study domain by accumulating the deviations from the mean and recording the absolute maxima and minima. This method is generally classified as a visual technique where changes in the slope of the linear trend are attributed to changes in mean values of the time series. Following the approach of Yousef, et al. [24], we rely on the peaks and troughs of the CUSUM plot, rather than slope changes, to indicate years where the most abrupt deviation from the of mean occurs. Then, the spatial distributions of the obtained change points (years) are mapped and examined.

Spatial Analysis
First, we focus on the comparison of the spatial organization of the three datasets plotted at their native resolutions. The time-averaged spatial distribution of TWS, precipitation, and soil moisture are obtained for winter (DJFM) and summer (JJAS) months during their common coverage period from 2002 to 2011. Figure 2 shows the winter and summer spatial distribution of TWS (a,b), precipitation (c,d), and soil moisture (e,f) derived over the 10-year study period.
Subplots 1a and 1b show the TWS anomalies varying in the range of −6.2 to 3.2 cm/month and −6.8 to −0.3 cm/month for the winter and summer periods, respectively. For the winter period, the largest decreasing anomalies (−5 to −6.2 cm/month) are recorded in north-western Saudi Arabia and the Empty Quarter (Rub' al Khali) region, whereas the largest increasing anomalies (2.5 to 3.2 cm/month) are recorded along the western coast, particularly around Jeddah, and the southern coasts of Oman and Yemen. For the summer period, decreasing anomalies are recorded over the entire AP, with the largest declines in the northern quadrant near Riyadh (−5 to −6.8 cm/month), and the least declines near Oman (−1 to 0 cm/month).
During the winter period, subplot 1c shows the TRMM precipitation retrievals reaching over 100 mm in the northern parts of Saudi Arabia, while negligible amounts are observed over the southern parts of the AP. Conversely, for the summer period, subplot 1d indicates negligible precipitation over the central and northern parts of the AP, while localized and intense events recorded up to 117 mm/month over the Yemen coastal highlands, as well as light rainfall (20 to 40 mm/month) over inland Yemen and the eastern coast of Oman.
For both winter and summer periods, subplots 1e and 1f show the highest readings of soil moisture values ranging from 0.3 to 0.45 g/cm 3 /month along the western coast, southern coast of Yemen, and its highlands, as well as the eastern coast of Oman and the Hajar Mountains on the border with the UAE. However, larger areal extents of soil wetness (0.2 to 0.3 g/cm 3 /month) are observed for the winter period over the northern AP and areas between Saudi Arabia, Oman and UAE.
The examined datasets show reasonable similarity in the spatial distribution of water storage on one hand and rainfall, and soil moisture on the other. Areas where soil moisture and precipitation values are the highest, in terms of both magnitude and extent, correspond to the significant positive TWS anomalies. The direct dependence of soil moisture on precipitation is less evident, especially over the northern part of the AP during summer, where the sustained soil wetness is observed in the AMSR-E product despite the absence of rainfall retrievals from TRMM.

GRACE Trends over the UAE and Comparisons to Well Data
Historically, the northern part of the UAE is known to be highly susceptible to the majority of rainfall events across the country. It is also considered to be a significant recharge area due to the coupled effect of minimal urbanization and highly permeable soils. Hence, any variation in rainfall is expected to directly reflect on TWS fluctuations. The decreasing trend of well-derived anomalies shows seasonal and intra-annual agreement with the GRACE trend over the UAE (see Figure 3). The sharp decrease in the GRACE anomalies to −5.1 cm/month in July 2002 GRACE is attributed to tapering artefacts during the initial 6-month calibration and validation period (April-September, 2002) which is excluded from the following analyses [49,50]. varying throughout the years. After 2010, both time series show a consistent decrease with low seasonal variability, failing to recover to positive anomaly magnitudes. Both time series resulted in significant negative trend coefficients at the 5% significance level. Based on the established agreement, GRACE anomalies are used to conduct the study over the regional scale of the AP.

Modified Mann-Kendall Test and Spatial Linear Trends
The linear trend coefficients are calculated for each dataset during the overlapping coverage (2002-2011) and mapped in Figure 4 for (a) GRACE, (b) TRMM, and (c) AMSRE. The Modified Mann-Kendall Test was applied to the derived linear trends to assess their significance at a 5% significance level. The TRMM daily rainfall trends were found insignificant at the 5% significance level, but held significant at the 10% level. GRACE monthly water storage anomaly trends were found significant across the entire AP. Increasing AMSR-E monthly soil moisture trends in the northern regions were significant, while the decreasing trends in the southern regions were insignificant at the 5%, but tested significant at the 10% level. The time series agree in terms of annual peaks and their consistent occurrence in winter time, namely, December and January. However, a slight shift is observed between the two time series with the GRACE anomalies consistently peaking approximately 2 to 4 months before the well-derived anomalies. This is primarily explained by the difference in scale between the GRACE grid size (~110 km × 110 km = 12,000 km 2 ) and the sub-scale area covered by the UAE well network (~1000 km 2 ) depicted in Figure 1c. The point-based well observations are sensitive to localized rainfall events over the network area, whereas GRACE detects total water column variations over its 12-fold larger grid area subject to larger scale rainfall systems.
Since GRACE detects changes in the total water column, i.e., subsurface storage + soil moisture + snow + canopy + river storage, the well-derived anomalies must be exclusively compared to the subsurface component. However, the snow, canopy and river storage components are not relevant to the study area. Moreover, the high operational frequency of ASMR-E (10.7 to 36.5 GHz) soil moisture acquisitions limits their penetration depth to the top skin layer of the soil which is subject to record-high evaporation rates [35]. Additionally, the high soil porosity in the region implies that the surface moisture either rapidly evaporates or infiltrates to deeper soil layers, which AMSR-E frequencies are not sensitive to [3]. This is evident in Figure 2 where the AMSR-E soil moisture distribution is largely similar between the wet and dry seasons. The magnitude of the peaks has been varying throughout the years. After 2010, both time series show a consistent decrease with low seasonal variability, failing to recover to positive anomaly magnitudes. Both time series resulted in significant negative trend coefficients at the 5% significance level. Based on the established agreement, GRACE anomalies are used to conduct the study over the regional scale of the AP.

Modified Mann-Kendall Test and Spatial Linear Trends
The linear trend coefficients are calculated for each dataset during the overlapping coverage (2002-2011) and mapped in Figure 4 for (a) GRACE, (b) TRMM, and (c) AM-SRE. The Modified Mann-Kendall Test was applied to the derived linear trends to assess their significance at a 5% significance level. The TRMM daily rainfall trends were found insignificant at the 5% significance level, but held significant at the 10% level. GRACE monthly water storage anomaly trends were found significant across the entire AP. Increasing AMSR-E monthly soil moisture trends in the northern regions were significant, while the decreasing trends in the southern regions were insignificant at the 5%, but tested significant at the 10% level. Remote Sens. 2021, 13, x FOR PEER REVIEW 11 of 18 Overall negative trend coefficients ranging from −8.6 × 10 −2 to −0.5 × 10 −2 cm/month are obtained across the entire study domain for GRACE, while both positive and negative trend coefficients ranging from −2.3 × 10 −4 to 3.2 × 10 −4 mm/day and, −1.5 × 10 −4 to 2.7 × 10 −4 g/cm 3 /month are determined for TRMM and AMSR-E, respectively. Overall negative trend coefficients ranging from −8.6 × 10 −2 to −0.5 × 10 −2 cm/month are obtained across the entire study domain for GRACE, while both positive and negative trend coefficients ranging from −2.3 × 10 −4 to 3.2 × 10 −4 mm/day and, −1.5 × 10 −4 to 2.7 × 10 −4 g/cm 3 /month are determined for TRMM and AMSR-E, respectively.
The largest negative trend coefficients for GRACE (−8.6 × 10 −2 cm/month) are located in the northern part of the AP. Figure 4b,c show decreasing trends for TRMM and ASMR-E in the southern part of the AP and increasing trends in the northern part. There is a noticeable agreement between TRMM rainfall and AMSR-E soil moisture trends and their spatial distribution, which is not captured in the analysis of Figure 2. Although GRACE data show a negative trend across the entire study domain, the most significant decline is recorded in the northern regions, where, conversely, rainfall and soil moisture experience increasing trends. Further analysis involved dividing the time series to pre-and post-change point years, and assessing whether the slope directions showed different results than those obtained from the entire coverage period. This analysis was carried out for TRMM alone since GRACE had an insufficient temporal coverage for time series division. A sample size larger than 199 entries is recommended for the application of the CUSUM method [51,52]. A major change in precipitation variability in the region was marked by the occurrence of change points. Increasing trends are obtained for both the pre-and post-change point time series across the entire study domain. Recall that this finding is not consistent with the outcome of the analysis for the entire time series with decreasing trends over most of the southern part of the AP. This indicates that the change in the time series around the detected break point is significant and is associated with a change in the overall precipitation regime in the region which is discussed below. Remote Sens. 2021, 13, x FOR PEER REVIEW 13 of 18 Further analysis involved dividing the time series to pre-and post-change point years, and assessing whether the slope directions showed different results than those obtained from the entire coverage period. This analysis was carried out for TRMM alone since GRACE had an insufficient temporal coverage for time series division. A sample size larger than 199 entries is recommended for the application of the CUSUM method

Spatial Analysis
The western coast of the AP and northern parts of Saudi Arabia record the highest values of the three key variables, particularly for the winter season, while higher values of TWS and precipitation are observed in the southern mountainous regions during the summer period. These results are related to the drastically different synoptic dynamics characterizing the northern and southern part of the AP. In fact, the northern and central regions of the AP receive the majority of precipitation during the winter season, in connection with the intrusion of extra-tropical cyclones from the Mediterranean, while the southern regime is mostly regulated by the summer tropical storm activity of the Arabian Sea and orographic processes [5,24,53].
Soil moisture monitoring is a challenging task, particularly in arid and hyper-arid environments. The persistently high soil moisture retrievals along the coastal areas are attributed the mixed sea-land pixels, while the higher soil moisture values in the north are attributed to the irrigation of agricultural establishments along with winter precipitation. Moreover, the high summer soil moisture values over the Hajar Mountains along the UAE border are the result of isolated orographic rainfall events of limited geographic extent triggered by local microscale weather processes [54]. Such local and short-lived rainfall events, although noticeable in the AMSR-E soil moisture product, are not captured by the TRMM product due to its coarse spatiotemporal resolution (the 0.25 • daily product was used here), which does not necessarily coincide with their occurrence. It is also worth noting that the retrieval of soil moisture using passive microwave observation is less reliable over mountains and complex terrain [35]. This can undermine the quality of the product over the highlands in the southwestern region. Another source of uncertainty in the soil moisture retrieval in the region is attributed to the discrepancy between the passive microwave signal in terms of brightness temperature and the thermal temperature required for the retrieval. Both temperatures have different penetration depths in desert land cover. The thermal temperature is sensitive to the top few millimeters of soil whereas the microwave brightness temperature penetrates deeper in the soil layer depending on the soil texture [55,56]. The lack of in situ soil moisture observations and the challenges in carrying out such observations in the region hamper the verification and the enhancement of the local retrievals [12,57,58].
Despite negligible summer precipitation over the northern AP, the soil wetness is sustained, which is explained by the increasing agricultural activity in those areas [19]. More importantly, these areas correspond to the largest TWS depletion anomalies, which indicates the excessive reliance on groundwater extraction for irrigation. Infiltration of the northern winter precipitation marginally compensates for the groundwater extraction, as suggested by the stable anomalies (close to null). At this stage of the analysis, the results corroborate findings from previous studies regarding the control of anthropogenic factors, namely groundwater extraction, on water storage depletion over the northern AP [18]. However, the impact of climatic changes, particularly precipitation, on TWS variations across the domain remains less evident and is the subject of the following sections.

Linear Trend Analysis
The GRACE TWS depletion rates are found to be the highest (−6.9 × 10 −2 to −8.6 × 10 −2 cm/month) over northern Saudi Arabia and Kuwait, including, as expected, the aforementioned agricultural demand, which relies on excessive groundwater extraction. These depletion rates override the increasing trends of rainfall and soil moisture in that region in line with the findings of Sultan, et al. [18]. In fact, large areas of desert lands around central Saudi Arabia have been developed into highly productive farms with an actual agricultural water demand of 29.82 billion m 3 or 94% of the total water consumption in 1994 [59].
The decrease of TWS anomalies in the southern part-reaching −2 × 10 −2 to −0.5 × 10 −2 cm/month over Yemen and Oman-was not as significant as in the case of the northern part of the AP, despite the southern decreasing trends of rainfall and soil moisture. The milder climatology and lower evapotranspiration rates over this region permit more recharge and, consequently, less groundwater depletion rates. Additionally, groundwater migration paths over the AP were studied by Alsharhan, et al. [60] who reported macro-flows from the northeast of Saudi Arabia towards the southeast. This may favor a time-lagged positive feedback between the increasing northern trends of rainfall (and soil moisture) and the southern TWS variability [3].
Oil production in the area, which requires pumping considerable volumes from underground reservoirs, may also be reflected in the storage variability. However, this storage loss is not considered in this work as most oil production activities are carried out in offshore locations. Furthermore, the seasonality of TWS variability as indicated by GRACE values (Figure 3) suggests that oil production does not dominate storage variability.

Change Point Monitoring
Considering the entire time series for assessing the trend direction and magnitude may be erroneous, particularly in the presence of change points in the time series marking a significant change in trend. Recall that it was noticed, in the case of TRMM data, which post-and pre-change point trends were increasing, while the overall time series showed a decreasing trend. Based on the comparison of the spatial distribution of change points between TRMM and GRACE, it is evident that there is a consistent pattern over central Saudi Arabia and Yemen. This implies that TWS variability at the macroscale is highly susceptible to changes in precipitation, and not solely controlled by anthropogenic extraction. Increasing trends of rainfall (pre-and post-2005 time series) were obtained over the southern part of the AP as opposed to the decreasing trends of the undivided time series. The increasing rainfall trend changes over the UAE agree with the increasing trends reported by Ouarda et al. [23] from their analysis based exclusively on rain gauge observations. They linked changes in the UAE rainfall regime to the variability in the equatorial Pacific sea surface temperatures-a teleconnection triggered by shifts of the upper level stream towards the Equator during the positive El Niño phase [61], ultimately impacting moisture in the region [62].
The difference in overpass time with AMSR-E is augmented by the high potential evapotranspiration in the region, which causes most of the rain to evaporate from the skin soil layer following a rain event. Additionally, ASMR-E uses high frequency, not L-band, which dictates low penetration depth and more sensitivity to the skin soil layer subject to the high evaporation. Nevertheless, soil moisture is an instantaneous short lasting response because of rapid drainage and evaporation. El-Nesr, et al. [63] explored changes in the reference evapotranspiration using weather station data across the AP from 1980 to 2008 to identify trend directions as an indicator to climate variability in the region. The annual daily evapotranspiration average was shown to increase from about 9.6 mm day −1 in 1980 to 10.5 mm day −1 in 2008. The increasing evapotranspiration trends explains the inconclusive change point statistics of the AMSR-E product, as most of the rainfall is being recaptured by the atmosphere rather than retained as soil moisture and partitioned to the subsurface. Future work will involve using multi-satellite observations to generate longer time series for a more reliable analysis of the trends. Water storage variability could be extended through the use of surface and groundwater models as water storage observation are not available earlier. This is particularly challenging in the absence of reliable and comprehensive geological surveys in the region which makes the use of remote sensing as proposed in this study an attractive alternative.

Conclusions
The AP is known to be under severe water stress with only 200 m 3 /cp/year of renewable water which is alarmingly below the conventional water stress threshold of 500 m 3 /cp/year [64]. Hence, understanding the macro processes by large-scale monitoring of key hydrological variables, namely, rainfall, water storage, and soil moisture is essential. The time series of these parameters are obtained from remotely sensed datasets, namely, TRMM, GRACE, and AMSR-E. Trend analysis reveals an overall depletion of TWS over the entire AP, especially, the northern part while rainfall and soil moisture trends were decreasing over the southern part. Alternatively, an overall increasing rainfall trend was observed after change point division of the TRMM time series. Change point detection revealed significant spatial dependencies between the TRMM and GRACE products, especially around Yemen during 2010, indicating hotspots of climate change impacts on TWS variations that should be accounted for alongside anthropogenic factors.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.