Assessing the large-scale plant-water relations using remote sensing products in the humid subtropical Pearl River Basin in south China

Vegetation interact closely with water resources. Conventional studies of plant-water relations at the field scale are fundamental for understanding the mechanisms of how plants alter and adapt to environmental changes, while large-scale 10 studies can be more practical for regional land use and water management towards mitigating climate change impacts. In this study, we investigated the changes in total water storage (TWS), aridity index (AI) and vegetation greenness, productivity and their interactions in the Pearl River Basin since 2002. Results show overall increase of TWS especially in the middle reaches where vegetation greenness and productivity also increased. This region dominated by croplands was identified as the hotspot for changes and interactions between water and vegetation in the basin. Vegetation was more strongly affected 15 by TWS than precipitation (P) at both the annual and monthly scales. Further examination showed that the influence of P on vegetation in wet years was stronger than dry years, while the TWS impact was stronger in dry years than wet years; moreover, greenness responded faster and productivity slower to dryness changes in dry years than wet years. The lag effects resulted in nonlinearity between water and vegetation. This study implies that vegetation in the basin uses rainwater prior to water storage until it gets dry, and the degree of water restriction on vegetation was higher than that of water consumption by 20 vegetation even in this rain-abundant region.


Data sources, pre-processing and analysis
To assess the plant-water relations at a large spatial scale, we obtained hydroclimate and vegetation data from different sources. Precipitation (P) and potential evapotranspiration (ETp) were obtained from Global Land Data Assimilation System (GLDAS) (Rodell et al., 2004); aridity index (AI) was then calculated as the ratio of ETp to P to represent the dryness 100 condition. Total water storage (TWS) change is inferred by the mass change detected by GRACE satellites (Tapley et al., 2004). GRACE data can be accessed from the Jet Propulsion Laboratory (JPL), the Center for Space Research (CSR), and the German Research Centre for Geosciences. Previous studies have shown that the ensemble mean of different products is effective in reducing the noise in the gravity field solutions (Long et al., 2017;Sakumura et al., 2014). Here we used total water storage anomaly (TWSA) data from the JPL and CSR with 'mascons' solution (release 6) at a resolution of 0.5° and 105 monthly. Cubic spline interpolation was applied to estimate the missing monthly data for the GRACEJPL and GRACECSR products during 04/2002-03/2015 that cover 13 hydrological years.
Vegetation data in this study include Normalized Difference Vegetation Index (NDVI) and Gross Primary Production (GPP).
NDVI was obtained from the GIMMS project during 04/2002-03/2015 and resampled to 0.5° using the nearest neighbor method and averaged to monthly to match the spatiotemporal resolution of GRACE and GLDAS data. Monthly GPP was 110 obtained from the Numerical Terradynamic Simulation Group in the University of Montana and rescaled to 0.5° (Running et al., 2004). Information of data sources, resolution and time span for all variables related to this study is listed in Table 1. To compare with GRACE data, anomalies of P, AI, NDVI, and GPP data were calculated by subtracting the means over the same baseline period of GRACE data (i.e. 01/2004-12/2009).
To investigate the changes in hydroclimate and vegetation, we carried out trend analysis using the Mann-Kendall (MK) test 115 method both in space and in time. The MK test does not require normality of time series and is less sensitive to outliers and missing values (Pal and Al-Tabbaa, 2009). This non-parametric test method has been used in many studies to detect changing hydrological regimes (Déry and Wood, 2005;Zhang et al., 2009). Interplay between hydroclimate and vegetation was quantified by linear regression; the Pearson correlation coefficient (r) and coefficient of determination (R 2 ) were taken as a measure for assessment of the linkages between different variables. Furthermore, a lag effect analysis was carried out to 120 determine the temporal dependency between variables where the linear relationship was not obvious.
Since the interactions between hydroclimate and vegetation can be different under dry and wet conditions, we hereby selected dry and wet years according to the anomalies of TWS, NDVI and AI under the criteria that dry conditions correspond to negative anomaly values of TWS and NDVI in addition to high anomaly of AI. Then the relationships between hydroclimate dynamics and vegetation greenness and productivity were specifically analysed. Uncertainties of the data used 125 were estimated by the standard error of each variable at the monthly and annual scales.

Changes in water storage and dryness
Comparison of the GRACE data from JPL and CSR shows that mean annual TWSA from GRACEJPL was overall greater than that from GRACECSR ( Fig. 2a-b). Both products showed clear zonal characteristics similar to the average of the two (Fig.  130 2c) that TWSA was generally higher in the middle-to-east areas than the rest of the basin especially the west upland, which infers a generally wetting condition in comparison to the baseline period. The trends of annual TWSA (Fig. 2d) showed that over the 13 hydrological years the TWS in most of the basin has increased at a rate below 10 mm yr -1 with 46% in the range of 5.0-10.0 mm yr -1 . Areas with low changing rate were mainly located in the west upland where the predominant land cover is grassland. Like the distribution of TWSA, water storage increase rate was also higher in the middle-to-east areas, where 135 overlap partly with croplands, than the rest of the basin.
Temporally, the basin has been getting wetter in general from 2002 (Fig. 2e). The TWSA has increased over the 13 years by 6.8 mm yr -1 inferred by GARCEJPL and 4.6 mm yr -1 by GRACECSR, with an average of 5.9 mm yr -1 . Noticeably, there were three shifts in the drying and wetting tendencies over the study period, i.e. the shift from drying between 2002 and 2005 to wetting between 2005 and 2008, followed by the shift back to drying between 2008 and 2011, and finally the shift to wetting 140 after 2011. In the following sections, only the mean TWSA from GRACEJPL and GRACECSR was used for analysis. Figure 3 shows the aridity index (AI) characterizing the spatial and temporal patterns of dryness. Majority of the basin has a semi-humid climate (AI=1.0~1.5); the west upland was clearly drier than the rest of the basin which is associated with precipitation patterns. Although dryness condition has not changed significantly over the 13 years with an overall positive trend in space (0.004±0.012) and time (0.007), it has some interesting characteristics such as the drying tendencies primarily 145 located in the southern areas, and the alternate periodical wetting and drying episodes temporally like TWSA. Areas with low TWS change rates coincided with drying climate represented by aridity index.

Changes in vegetation greenness and productivity
Spatial NDVI distributions (Fig. 4) were highly related to vegetation cover types that the high NDVI values coincided with forest covers and low values corresponded to impervious surfaces, grasslands and croplands. It clearly reflects the impacts of 150 urbanization on surface greenness particularly near the basin outlets in the southeast. Over the 13 years NDVI has not shown significant changes across the basin, since the majority (~70%) had a MK test p>0.05 at the pixel scale. The areas with significant changes were concentrated in the central south of the basin where croplands are predominant. This infers an intensification of crop farming activities over these areas.
Temporally, NDVI has an overall insignificant increase trend over the 13 years at an annual rate of 0.004 (p=0.56) with 155 interannual fluctuations. It is noticeable that the periodical shifts in the NDVI trends were almost identical to TWSA in Fig. 2e. This reflects a tight bound between the vegetation greenness and water availability in this rain-abundant region.
Interestingly, in 2004 when water storage continued to decrease following the previous years, NDVI did not show a https://doi.org/10.5194/hess-2020-242 Preprint. Discussion started: 29 June 2020 c Author(s) 2020. CC BY 4.0 License. continuity of decreasing but increased instead, implying a vegetation resilience and recovery after previous dry period. The recovery coincided with a slight decrease in aridity index, hence, vegetation did not respond solely to water availability but 160 also to atmospheric demand.
In addition to NDVI, GPP was also analysed for the basin (Fig. 5). It is not surprising to observe that GPP was highly responsive to NDVI such that areas with low NDVI also had low GPP (e.g., the central agricultural region and upland grassland). GPP anomaly also showed positive high values in the central south areas dominated by croplands coincident with NDVI anomaly, indicating an increased agricultural production induced by intensified agricultural activities in this region. It 165 should be noted that most of the trends were not statistically significant. Over the entire basin, annual GPP showed almost the same periodical decreasing and increasing trends as NDVI and TWSA, except that the third turning point occurred in 2010 rather than 2011. Linear regression gave a coefficient of determination R 2 =0.59 (p=0.002) between annual TWSA and NDVI, higher than that between TWSA and GPP (R 2 =0.12, p=0.257), which may imply a more direct and stronger impact of water stress on vegetation greenness than productivity at an annual scale. 170

Interactions between hydroclimate and vegetation
Combining Fig. 2-5, we found that climate condition, water storage and vegetation dynamics are tightly interlinked.
Spatially, precipitation, water storage and dryness affected vegetation in a similar way compared to temporal characteristics, i.e. the influence of TWS was relatively stronger than P and AI. The hotspots of the interactions were found in the middle areas, and dryness more negatively affected greenness than productivity in these areas (Fig. 7). These analyses indicate that atmospheric stress and water stress imposed more direct and stronger impact on vegetation greenness than productivity on a 180 yearly basis, and water constraint on vegetation was stronger than that of dryness.
At the monthly scale, however, the linear responses of GPP to P and TWS were stronger than the linear responses of NDVI to P and TWS ( Fig. 8a-b). The response of both NDVI and GPP to P was more nonlinear than to TWS, and the sensitivity of NDVI and GPP to TWS was stronger than to P indicated by the regression slopes, implying a stronger link between water storage and vegetation growth. Meanwhile, increase in dryness resulted in nonlinear decreases in NDVI and GPP (Fig. 8c). 185 The relationships show that although precipitation is the main water input to the terrestrial hydrological cycle, it is how much water is stored in the soils that determines vegetation greenness and biomass production at a shorter time scale than annual.
Nonlinear plant-water relationships can be explained by the lag effect that monthly changes of NDVI and GPP fell behind the changes of P and TWS to varying degrees (Fig. 9). This means that the water restriction on vegetation outweighed the water consumption by vegetation. Vegetation response to hydroclimate changes is expected to differ in dry and wet years. 190 Here, we assumed that the annual anomalies of TWS<0, NDVI<0 and AI>0 corresponded to dry conditions, and hence https://doi.org/10.5194/hess-2020-242 Preprint. Discussion started: 29 June 2020 c Author(s) 2020. CC BY 4.0 License. defined 2003, 2007and 2011as dry years and 2002, 2006, 2010 as relatively wet years. There was evidence of drought occurrences in these dry years (Lin et al., 2017;Wang et al., 2014b). It can be seen that the dry and wet years were mainly differentiated by the rainfall data in summer months July and August affecting water storage and dryness. The range of NDVI and GPP was 32.8% and 8.4% higher on average in dry years than wet years, mainly 195 attributable to the difference in the non-growing seasons from October to March. Both the minimum and maximum NDVI were lower in dry years than in wet years, particularly, the minimum NDVI in dry years was 87.4% lower than that in the wet years, compared to 8.6% lower for the maximum. GPP was similar in dry and wet years, with only 9.8% and 6.9% lower in dry years for minimum and maximum values, respectively. This implies that vegetation greenness is more sensitive to any changes in hydroclimate than productivity. Moreover, GPP in non-growing seasons in dry years was relatively higher than 200 that in wet years reflecting a positive effect of stress on biomass accumulation. Figure 10 gives the R 2 from linear regression between different variables considering phase shift for lag analysis. It shows NDVI varied strongest with P, TWSA and AI in the previous 3, 1 and 3 months, respectively when considering all data during 2002-2014. In comparison, a shorter lag time of GPP to P, TWSA, and AI was detected (2, 0, 1 month, respectively).
Comparison of the lag time in dry and wet years shows that the influence of P on vegetation was more prominent in wet 205 years than in dry years, while TWS influence was greater in dry years than wet years. Moreover, NDVI responded faster to dryness change in dry years (2 months) than wet years (3 months), and GPP responded slower to dryness change in dry years (1 month) than wet years (0 month). This may indicate that dryness can stimulate biomass production to some degree. In addition, GPP varied synchronously with TWS showing a high dependency on water storage despite the dryness conditions.

Uncertainties in the datasets and results
Data availability is one of the greatest obstacles for large-scale and long-term ecohydrological studies. Remote sensing products are thus useful to characterize ecohydrological changes in a large sparsely monitored basin. In this study, we used remote sensing and assimilated data of water storage, vegetation status and precipitation to assess their relationships.
GLDAS uses meteorological forcing data merged from multiple sources including ground and satellite observations, and 215 GLDAS precipitation proves to be highly consistent with observations in China (Mo et al., 2016;Wang et al., 2016). Here we also compared the GLDAS P with the measured P in the pixels where stations are available (Fig. 11). Overall, P from GLDAS agreed well with observations with R 2 ranging from 0.69 to 0.89 (±0.05) spatially, while on average the monthly P from GLDAS underestimated observations by ~10% over all valid pixels. The comparison provides some confidence in applying the GLDAS products for long-term hydrological trend analysis, though discrepancies exist in the absolute values. across China and concluded that MODIS GPP was more reliable over grassland, cropland and mixed forestland than the other datasets. These land cover types happen to be the predominant ones in the Pearl River Basin. Zhang et al. (2017b) and Yuan et al. (2015) also compared various GPP datasets globally and regionally, and inconsistencies existed in these comparisons that could stem from the way each algorithm parameterizing atmospheric and water stress and difference in the vegetation index data (Yuan et al., 2015). Despite the dispute of data accuracy, MODIS GPP seems more frequently used 240 due to its moderate spatiotemporal resolution and data coverage.
Inspired by the studies of TWS change using GRACE satellite data with different processing algorithms (Long et al., 2017;Sakumura et al., 2014), it may be more accurate and informative by using the average values from as many available datasets for the targeted ecohydrological variables as possible, i.e. the ensemble means, than using a single dataset. This is worth further investigation which could enhance the studies in many ungauged basins for critical hydrological assessments. 245

Hotspot for hydroclimate and vegetation changes
NDVI and GPP shared the same spatial patterns and high GPP corresponded to high NDVI in the forested areas. Low values existed in the west upland with grass cover and the central south areas of croplands. Over the 13 years NDVI and GPP showed insignificant changes with large interannual variabilities. Unlike the north China where vegetation cover is deeply affected and largely recovered through decades of ecological restoration projects (Chen et al., 2019;Feng et al., 2005), 250 vegetation cover especially the forest cover which occupies most of the PRB almost remained constant from early 2000s . Even so, we identified the areas with significant increase in NDVI and GPP in the central south region of the basin where croplands dominate. Therefore, considering that the precipitation gradually decreases from southeast coastal area toward northwest outback of the basin, changes of TWS, NDVI and GPP jointly imply that the water storage increase in this hotspot region has resulted in the intensification of agricultural activities and boosted the food production since the early 2000s. It is for the first time in studies to reveal such phenomenon and can be meaningful for the food-water nexus studies in this region, and informative for a possible shift of China's main food production from the north to the south in the context of water richness in the south and shortages in the north (Kuang et al., 2015).

Causal roles of water and dryness in vegetation changes
The overall TWS increase is promising for the managers and users of water resources in the PRB, however, the strong 260 correlation with precipitation seasonality restrained the available water in the relatively dry periods. In fact, previous studies have reported the contribution and restriction of P to TWS. For instance, Chen et al. (2017) revealed the liability of P to TWS (r=0.78) in the PRB. Mo et al. (2016) found TWS more strongly explained (60%) by annual P in river basins in south China than in north China. In this sense, storage shortage in dry periods subject to seasonal reduction of precipitation would hamper vegetation growth. Analysis in this study shows that NDVI was highly correlated with TWS and P at the annual 265 scale (Fig. 6), consistent with previous studies in the PRB and other areas (Guan et al., 2015;Zhaos et al., 2016;Zhu et al., 2018). Whilst at the monthly scale NDVI was still strongly influenced by TWS but not by P, in comparison to the strong response of GPP to both P and TWS. The weakened linear influence of P on NDVI at the monthly scale, found also by others such as Bai et al. (2019) andA et al. (2017), can be explained by the lag effect that NDVI lagged by 3 and 1 months after P and TWS, respectively. In comparison, the lag time between GPP, P and TWS was 1 month shorter than NDVI 270 versus P and TWS. In addition, comparison of the plant-water relations in dry and wet years showed a slower response of GPP to aridity index in dry years than wet years (Fig. 10b-c), which may imply that a certain degree of drying can stimulate biomass accumulation. This phenomenon is similar to the principle of regulated irrigation in agriculture to increase water use efficiency under a certain degree of water stress (Chai et al., 2016), and also revealed by other studies . This dryness effect on ecosystem productivity cannot be detected in the annual scale assessment (Brookshire and 275 Weaver, 2015;Yao et al., 2020). These results indicate firstly that pre-growing season hydroclimate conditions play a key role in the follow-on vegetation growth and production , and secondly that water limits vegetation even in this subtropical rain-abundant region instead of water shortage resulted from vegetation establishment.
Anomalies of TWS, aridity index and NDVI together well defined the occurrences of drought in the basin that are identical to other studies using P, TWS alone or other drought indices (Wang et al., 2014b;Zhang et al., 2018). The drying episodes 280 confined the vegetation greenness and production (Lin et al., 2017). Liu et al. (2014b) reported that China's national total annual net ecosystem productivity exhibited declines during 2000-2011, mainly due to the reduction in GPP caused by extensive drought. Although drought is generally associated with declines in vegetation greenness and productivity due to water and heat stresses (Eamus et al., 2013), the magnitude of vegetation reduction, determined by ecosystem sensitivity to drought, can vary dramatically across plant communities. While Zhang et al. (2017a) detected insensitivity of vegetation to 285 droughts in humid south China including the lower reach of PRB, this study observed that NDVI experienced a recovery in 2004 after drought in the previous year, which may be a result of irrigation during drought in the agricultural regions since forests are more resilient to droughts (DeSoto et al., 2020;Fang and Zhang, 2019). Future climate projections predict https://doi.org/10.5194/hess-2020-242 Preprint. Discussion started: 29 June 2020 c Author(s) 2020. CC BY 4.0 License.
increases in temperature and insignificant changes in precipitation in the basin which would trigger more heatwave induced flash droughts . To mitigate the impacts on both water resources and ecosystems, proper plans should be 290 made such as conversion of the low resilient ecosystems to forests (Fang and Zhang, 2019) and improvement of biodiversity in ecosystems (Isbell et al., 2015;Oliver et al., 2015), in addition to engineering regulations like reservoir operations (Lin et al., 2017).

Conclusions
Plant-water relations over the Pearl River Basin were examined using remote sensing products during the hydrological years 295 of 2002-2014. Results show that water storage has increased across the entire basin at an average rate of 5.9 mm yr -1 .
Vegetation greenness and productivity has also shown some changes but not overall significant. Spatial characterization reveals that the central south areas of the basin dominated by croplands are the hotspots for the changes of and interactions between hydroclimate and vegetation. This implies an increase in food production induced by intensification of agricultural activities in these areas. Lag effect analysis at the monthly scale reflects that even in this rain-abundant subtropical basin the 300 water restriction on vegetation precedes the water consumption by vegetation. Furthermore, comparison of the plant-water relations in dry and wet years showed a stronger influence of precipitation and a weaker influence of water storage on vegetation in wet years than dry years. A slower response of vegetation productivity to aridity index in dry years than wet years was identified which may indicate a stimulating role of a certain degree of dryness on vegetation production. This study reveals the changes and interplay between plant and water using readily available remote sensing and assimilated data, 305 and has implications for proper measures regarding land use alterations to mitigating frequent drought impacts on water resources and ecosystems under a warming climate.

Data availability
The original data in the study are available from the links given in Table 1; the processed data can be obtained from the corresponding author on request.

Figure 7.
Pearson correlation coefficient between annual anomalies of (a-c) precipitation, total water storage, aridity index 565 and NDVI; and (d-f) precipitation, total water storage, aridity index and GPP.