Future Changes in High and Low Flows under the Impacts of Climate and Land Use Changes in the Jiulong River Basin of Southeast China

: Climate change and human activities have profoundly affected the world with extreme precipitation, heat waves, water scarcity, frequent ﬂoods and intense droughts. It is acknowledged that climate change will persist and perhaps intensify in the future, and it is thus meaningful to explore the quantitative impacts of these changes on hydrological regimes. The Jiulong River basin serves as an important watershed on the southeast coast of China. However, future hydrological changes under the combined impacts of climate change and land use change have been barely investigated. In this study, the climate outputs from ﬁve general circulation models (GCMs) under the Coupled Model Intercomparison Project Phase 6 (CMIP6) were corrected and spatially downscaled by a statistical downscaling method combining quantile mapping and machine learning. The future high-resolution land use maps were projected by the CA–Markov model with land use changes from the Land-Use Harmonization 2 (LUH2) as constraints. The future dynamic vegetation process was projected by the Biome-GBC model, and then, the future hydrological process under four representative concentration pathways and shared socioeconomic pathways (RCP–SSP) combined scenarios was simulated by a distributed hydrological model. Based on the copula method, the ﬂood frequency and corresponding return periods were derived. The results demonstrated that future precipitation and air temperature would continue to rise, and future land use changes would have different developing pathways determined by the designs in various SSP–RCPs. Under the combined impacts of climate and land use change, the total available water resources will increase due to increasing precipitation, and the high ﬂow and low ﬂow will both increase at three stations under the four SSP–RCPs. The annual 1-day maximum discharge is projected to increase by 67–133% in the last decade of the 21st century, and the annual 7-day minimum discharge is projected to increase by 19–39%. The ﬂood frequency analysis showed that the Jiulong River basin would face more frequent ﬂoods in the future. By the end of the 21st century, the station-average frequency of a historical 100-year ﬂood will increase by 122% under the most optimistic scenario (SSP126) and increase by 213% under the scenario of greatest regional rivalry (SSP370). We demonstrated that climate change would be the major cause for the increase in future high ﬂows and that land use change would dominate future changes in low ﬂows. Finally, we recommend integrated and sustainable water management systems to tackle future challenges in this coastal basin. variability, with the three stations increasing by an average of 144% under SSP126, 69% under SSP245, 87% under SSP370 and 64% under SSP585 in the last decade of the 21st century. The results demonstrate that the annual 1-day maximum, 3-day maximum and 7-day maximum discharges all reveal increasing trends in the future under all scenarios, which indicates that there may be more frequent floods with larger magnitudes in the future.


Introduction
With socioeconomic developments, climate change and human activities have profoundly affected the hydrologic cycle around the world [1]. These changes have led to changes in the severity, frequency and duration of extreme events [2], including extreme heat waves [3], intense droughts [4] and recurrent flooding [5]. These extreme events may have a negative impact on biodiversity [6], crop yields [7], human health [2] and on various other socioeconomic factors [8]. Coastal areas typically have dense populations, concentrated urban construction and developed economies that are highly vulnerable to extreme floods [9,10]. Thus, it is imperative to explore the potential changes in future flow regimes by considering both climate change and human activities under global warming and socioeconomic development to formulate adaptation and mitigation strategies, especially for coastal basins [11].
Numerous studies have detected possible future changes in runoff under the impacts of climate and land use changes, and then analyzed with great detail the potential changes in different hydrological regimes including high and low flows. Most of them adopted hydrological models with projected future climate forcing data from general circulation models (GCMs) under the Coupled Model Intercomparison Project (CMIP) [12][13][14][15][16]. Based on a global river routing model with outputs of 11 GCMs from CMIP5, Hirabayashi et al. [14] projected future changes in runoff on a global scale and concluded that a warmer climate would significantly increase future high flows in Southeast Asia, peninsular India, eastern Africa and the northern half of the Andes. Meanwhile, their scenario analysis revealed that global flood risks will increase with climate warming. Using a distributed semi-physically based rainfall-runoff model with climate projections from CMIP5, Alfieri et al. [17] focused on the high flow changes and projected future flood frequency in Europe under warmer future conditions and concluded that flood peaks with return periods above 100 years would, on average, double in frequency within 3 decades. Rottler et al. [18] detected and assessed future changes in runoff in the Rhine River under various climate scenarios in CMIP5 with full consideration of basin hydrological processes by a mesoscale hydrologic model and demonstrated increases in high flow due to future increases in antecedent precipitation and temperature. With a calibrated hydrological model and UKCIP98 climate change scenarios, Arnell [19] explored the impact of climate change on low flows in Britain, and concluded that the larger variability in future climate would lead to a decrease in low flows and the degree of the impact would be greater in upland basins.
In addition to climate change, a large number of studies have proven that land use change has a significant impact on flow regimes [1,[20][21][22]. Baker et al. [21] used the Soil and Water Assessment Tool (SWAT) to simulate the hydrological process of a small basin in Kenya under three different land use maps and explored the impacts of forest degradation and cropland expansion on the allocation of surface and groundwater recharge as well as ecological effects. Their results showed that deforestation would shift the allocation and have a negative ecological effect. With a process-based hydrological model, Li et al. [22] detected the effects of urbanization on various elements of the hydrological cycle and concluded that the increase in impervious surfaces and decrease in evapotranspiration caused by vegetation losses during rapid urbanization would lead to an increase in low flow, and eventually increase in water yield. Considering the impacts of land use change and climate change simultaneously, some studies projected future possible changes in hydrological processes at a watershed scale [1,[23][24][25][26]. Tavakoli et al. [27] explored the response of extreme flows under the impacts of climate change and urban growth with a hydrological model, and concluded that extreme low flows would decrease noticeably due to climate change and be less sensitive to urban development, while extreme high flows would increase greatly due to both climate change and urban development. Using a hydrological model and a land use simulation model, Chacuttrikul et al. [24] assessed the impacts of climate and land use changes on river discharge in a small basin in Thailand and showed that the increase in future precipitation would lead to an increase in both low and high flows, and that continuous conversion from forest to paddy fields would cause an increase in both high and low flows. In addition, they concluded that land use change dominated historical runoff changes, while future discharge changes would be mainly caused by climate change. Similarly, Yang et al. [1] adopted the SWAT model to project future hydrological processes and the CA-Markov model to simulate future land use changes in the Luanhe River basin in China and projected future possible changes in runoff under continuous afforestation and global warming. They illustrated that climate change would increase the high flows and lead to a slight decrease in low flow seasons, and future afforestation would cause a reduction in low flows. In addition, some of the increase in runoff due to climate change would be offset by reduced runoff caused by afforestation, but climate change would still play a dominant factor in the future. Although the impacts of climate change and land use change on hydrological processes have been evaluated by many previous studies, most of the land use change scenarios are designed independently and are not consistent with climate change scenarios designed by CMIP. This deficiency may lead to a large bias in accurately reflecting future possible changes in hydrological regimes.
The newly released CMIP6 has set a series of more reasonable scenarios than previous scenarios designed in CMIP5. The new scenarios are rooted in socioeconomic trajectories and are named representative concentration pathways (RCPs) and shared socioeconomic pathways (SSPs) combined scenarios [28]. Serving as a required input of GCMs in CMIP6, Land-Use Harmonization 2 (LUH2) provides a newly harmonized annual land use dataset from 850 to 2100 under future SSP-RCP scenarios [29]. To detect future changes in hydrological processes, a physically based hydrological model is necessary, and the geomorphologybased hydrological model (GBHM) developed by Yang et al. [30] can be used for future hydrological projections. The GBHM model is a fully distributed hydrological model and has been successfully applied in many regions [31][32][33][34][35] including the Jiulong river basin. However, there is a gap between the required inputs and the data provided by CMIP6 and LUH2, which have systematic bias and coarse spatial resolution. To bridge this gap, statistical downscaling methods can effectively correct systematic bias and downscale the spatial resolution and are widely used for future climate downscaling [36][37][38]. The CA-Markov model is a popular and robust land use simulation model and can potentially be adopted for spatial downscaling of land use maps [1,[39][40][41][42]. Combining the SSP-RCP scenario-based climate forcing from CMIP6 and land use maps from LUH2, it is possible to explore future water-related disaster changes considering greenhouse gas emissions and socioeconomic development simultaneously with the help of the GBHM and CA-Markov models.
In this study, we conducted climate downscaling using a statistical downscaling method and projected future high-resolution land use maps under various SSP-RCPs based on CA-Markov and LUH2. With future climate forcing and land use data, we projected future dynamic vegetation growth dynamics by a process-based ecosystem model and then simulated future hydrological responses in the Jiulong River basin. Based on these projections, we analyzed future changes in annual runoff and flow regimes which included high flows (annual 1-day maximum, 3-day maximum and 7-day maximum discharges) and low flows (annual 7-day minimum, 30-day minimum and 90-day minimum discharges), and calculated flood frequency using a copula-based method. We strived to (1) illustrate possible changes in future climate and land use; (2) demonstrate future changes in high and low flows under the combined impacts of climate and land use changes; and (3) provide corresponding suggestions to mitigate and adapt to future water-related risks.

Study Area
The Jiulong River basin (see Figure 1) originates in the Boping Mountains and is a coastal watershed with a drainage area of 14,741 km 2 . The river length is 285 km, and is the second longest river in Fujian Province, China. The basin is located in the subtropical monsoon climate zone with an average annual runoff of 1010 mm. The annual average temperature ranges from 19.9 • C to 21.1 • C, and the annual average precipitation is 1400-1800 mm, with approximately 75% of the annual precipitation occurring from April to September. The basin mainly consists of two tributaries, the North River and the West River, which cover nearly 95% of the total drainage area of the Jiulong River basin. The Jiulong River basin is profoundly impacted

Data
(1) Historical observation data The observed meteorological data used in this study mainly contained daily air temperature (mean, maxima and minima), daily relative humidity, daily sunshine duration and daily average wind speed, which were provided by the China Meteorological Administration (CMA). The gauged meteorological data were spatially interpolated to 1-km grids by an angular distance weighting method with elevation corrections [43]. Precipitation data were obtained from the 0.25° × 0.25° China Gauge-based Daily Precipitation (CGDPA) Product [44] and were resampled to 1-km grids by nearest neighbor assignment

Data
(1) Historical observation data The observed meteorological data used in this study mainly contained daily air temperature (mean, maxima and minima), daily relative humidity, daily sunshine duration and daily average wind speed, which were provided by the China Meteorological Administration (CMA). The gauged meteorological data were spatially interpolated to 1-km grids by an angular distance weighting method with elevation corrections [43]. Precipitation data were obtained from the 0.25 • × 0.25 • China Gauge-based Daily Precipitation (CGDPA) Product [44] and were resampled to 1-km grids by nearest neighbor assignment in ArcGIS Software. Four historical land use maps in 1980,1990,2000 and 2010 were obtained from the Geographical Information Monitoring Cloud Platform (http://www.dsac.cn/DataProduct/Detail/20080422, accessed on 15 January 2022) at a 30-m resolution and were reclassed into 10 types: water body, urban area, bare soil, forest, irrigated cropland, upland, grassland, shrub, wetland, and ice. Daily streamflow observations at three gauges (i.e., Zhangping, Punan and Zhengdian) were obtained from the Hydrology and Water Resources Bureau. Historical LAI data were obtained from the 3rd generation of the vegetation leaf area index product (LAI3g) of the Global Inventory Modeling and Mapping Studies (GIMMS) [45]. The soil property and soil hydraulic parameter data were obtained from the China Soil Hydraulic Parameters Dataset by Duan et al. [46]. The digital elevation data were obtained from the shuttle radar topography mission (SRTM) [47] database with a resolution of 90 m.
(2) Climate outputs of GCMs and land use projections from LUH2 In this study, five newly released GCM outputs (i.e., ACCESS-ESM1-5, GFDL-ESM4, MIROC6, MPI-ESM1-2-LR and MRI-ESM2-0) under four SSP-RCP scenarios were selected from the latest CMIP6. SSP1-RCP2.6 is referred to as SSP126, SSP2-RCP4.5 is referred to as SSP245, SSP3-RCP7.0 is referred to as SSP370, and SSP5-RCP8.5 is referred to as SSP585. The GCMs selected in this study have all been applied in Asia for future climate change research in previous studies [48][49][50][51][52][53], and detailed information is shown in Table 1. The climate model outputs from 1985 to 2100 mainly contain the daily near-surface mean, maximum and minimum air temperature, daily precipitation, daily near-surface relative humidity and daily near-surface wind speed. Rooted in the socioeconomic trajectories, the newly proposed SSP-RCPs can be harmonized with the RCP previously developed in CMIP5 with shared policy assumptions and provide more plausible future scenarios [28]. The SSP1 scenario represents a green-fueled growth future; SSP2 represents a 'middle-ofthe-road' future; SSP3 assumes high inequality between countries; and SSP5 represents a fossil-fueled growth future [54]. Future land use data were obtained from the LUH2 dataset, which provided harmonized annual land use states and transitions with a spatial resolution of 0.25 • × 0.25 • under four SSP-RCPs from 850 to 2100. Detailed information on LUH2 can be found in the study by Hurtt et al. [29]. Compared with observations, the climate outputs from GCMs always contain systematic bias, and their coarse spatial resolution makes them unable to meet the high requirements of hydrological models. Therefore, it is necessary to perform bias correction and spatial downscaling on GCM outputs. In this study, we used a statistical downscaling method that combines quantile mapping (QM) and a machine learning method. The framework of this method ( Figure 2) is summarized as follows: (1) The first step is the initial bias correction and spatial downscaling by quantile mapping (QM) [60] for the historical period and detrended quantile mapping (DQM) [61] for the future period. This step is applied to all climate variables. By matching the cumulative distribution functions (CDFs) of GCM outputs with those of observations, QM has been widely used to correct the bias between GCM outputs and observations in the field of climate change [36,38]. For bias correction in the future period, DQM removes the trend of simulations firstly and then conducts quantile mapping, eventually restoring the trend. This method can transfer the future distribution to lie within the support of the historical distribution. A detailed description of QM and DQM can be found in the studies by Piani et al. [60] and Cannon et al. [61]. Since QM or DQM mainly corrects the daily simulations on a monthly basis and the daily variability of precipitation is relatively larger, QM/DQM alone cannot meet the high requirements of daily hydrological simulation. Therefore, other methods are still needed to improve the correction at a daily scale. details can be found in the IDIRISI manual [69]. In this study, we use the future land use transition area matrixes from the LUH2 dataset to replace the transition area matrixes generated by the Markov process and then run the CA-Markov model to project future high-resolution land use maps under four SSP-RCPs. The model inputs mainly include the historical land use maps (100 m × 100 m) in 2000 and 2010, and the model finally generates future 100 m × 100 m land use maps for each decade from 2020 to 2090. Step 1: monthly downscaling by quantile mapping Step 2: daily downscaling by machine learning Step 3: extreme value correction by TGDM The temporal and spatial dynamic distribution of vegetation is affected by both climate and land use and then constrains the hydrological process through stomatal conductance, leaf area, interception and transpiration [70][71][72][73][74]. Thus, the process of dynamic vegetation growth should be considered in hydrological simulations. As a one-dimensional process-based biogeochemical model, the Biome-BGC model [75] can robustly and accurately describe vegetation physiological processes, including photosynthesis, respiration, evapotranspiration and decomposition, in terrestrial ecosystems. The leaf area index (LAI) is estimated as a function of leaf carbon amount and is updated daily [76]. The details of Biome-BGC can be found in the studies by White et al., [77] and Thornton et al., [75]. In this study, the 4.2 version of Biome-BGC is mainly used to project future LAI with future climate forcing and land use projections under four SSP-RCPs. The eco-physiological parameters mainly include maximum stomatal conductance, fraction of leaf N in Rubisco, canopy average specific leaf area, canopy water interception coefficient, current growth proportion, year day for the start of new growth and year day for the end of litterfall, and are calibrated with remote sensing-based LAI from 1985 to 1989. (2) The second step is daily correction by machine learning. This study selected the GA-NARX [32,62] method to perform further correction for precipitation. As a class of recurrent neural networks, GA-NARX is regarded as a suitable method for nonlinear dynamic time series prediction, and a detailed description can be found in [32,62]. In this study, the predictors of GA-NARX include the QM/DQM-corrected daily rainfall, mean/maximum/minimum temperature, relative humidity and wind speed, and the target is the observed precipitation. Similar to other machine learning methods based purely on data, GA-NARX always underestimates high values [63] due to the sigmoidal activation function and insufficient training with limited high values in the time series [31,32,64]. This drawback can significantly affect the projection of future floods and needs to be overcome.
(3) The third step is extreme value correction by triple-gamma distribution mapping (TGDM). TGDM belongs to the categories of QM methods that fit CDFs by a triple-gamma distribution similar to the double-gamma distribution in Gao et al. [65] and can better capture the distribution of extremely high values. In a triple-gamma distribution, the observations and simulations are separately divided into three groups according to the specifically designed percentiles, including the 95th percentile and 99th percentile, as criteria for distinguishing normal, high and extremely high values in this study. The detailed calculation process is similar to the double-gamma distribution in Gao et al. [65] and can be found in the corresponding reference.

Future Land Use Projection by CA-Markov and Leaf Area Index Projection by Biome-BGC
The CA-Markov model [41] is a popular method for modeling the dynamics of spatiotemporal land use with high accuracy [66][67][68]. The Markov chain is a random process describing the transition of a system from one state to another and has the stable and memoryless property which means the probability of a random process depends only on its current state and time. A cellular automaton (CA) can independently change its state according to its previous state and that of immediate neighbors, which includes the spatial relationship in its outputs. Integrating the advantages of the Markov chain and CA, CA-Markov can reallocate land use during iterations until reaching the area totals predicted by the Markov chain. The CA-Markov model is run using IDIRISI software, and more details can be found in the IDIRISI manual [69]. In this study, we use the future land use transition area matrixes from the LUH2 dataset to replace the transition area matrixes generated by the Markov process and then run the CA-Markov model to project future high-resolution land use maps under four SSP-RCPs. The model inputs mainly include the historical land use maps (100 m × 100 m) in 2000 and 2010, and the model finally generates future 100 m × 100 m land use maps for each decade from 2020 to 2090.
The temporal and spatial dynamic distribution of vegetation is affected by both climate and land use and then constrains the hydrological process through stomatal conductance, leaf area, interception and transpiration [70][71][72][73][74]. Thus, the process of dynamic vegetation growth should be considered in hydrological simulations. As a one-dimensional processbased biogeochemical model, the Biome-BGC model [75] can robustly and accurately describe vegetation physiological processes, including photosynthesis, respiration, evapotranspiration and decomposition, in terrestrial ecosystems. The leaf area index (LAI) is estimated as a function of leaf carbon amount and is updated daily [76]. The details of Biome-BGC can be found in the studies by White et al. [77] and Thornton et al. [75]. In this study, the 4.2 version of Biome-BGC is mainly used to project future LAI with future climate forcing and land use projections under four SSP-RCPs. The eco-physiological parameters mainly include maximum stomatal conductance, fraction of leaf N in Rubisco, canopy average specific leaf area, canopy water interception coefficient, current growth proportion, year day for the start of new growth and year day for the end of litterfall, and are calibrated with remote sensing-based LAI from 1985 to 1989.

Hydrological Simulation by GBHM
As a fully distributed hydrological model, the geomorphology-based hydrological model (GBHM) proposed by Yang et al. [30] mainly consists of two parts: hillslope unit-based hydrological simulations and river routing along river networks. The hydrological simulation in the GBHM mainly includes canopy interception, evapotranspiration, surface flow, infiltration, unsaturated flow and groundwater flow, and river routing is based on the kinematic wave method. A detailed description of the GBHM can be found in the studies by Yang et al. [30] and Yang et al. [78]. In this study, the GBHM is used for natural hydrological simulation and projection in the Jiulong River basin with a spatial resolution of 1 km × 1 km. Since the model had been calibrated and validated with the Jiulong River basin in our previous study by Lu et al. [33], this study directly used the developed model to simulate the hydrological process during the historical period (1 January 1985-31 Decmber 2014) and the future period (1 January 2015-31 Decmber 2099). For the detailed GBHM development process in the Jiulong River basin, please refer to the previous study by Lu et al. [33].

Statistical Approaches for Analyzing the Changes in High and Low Flows
For the definitions of high flow and low flow, there are two common methods to describe them. One is to determine the high and low flows by the high and low percentiles based on the flow frequency distribution curve, and is widely used to validate the performance of hydrological models and can reflect the most common situations of high and low flows [32,79]. The other method of description is the annual maximum or minimum of specific moving window mean discharges, and is commonly applied to reflect flood or drought risks [27,80]. Since the changes in future hydrological hazards are focused more in this study, we selected the second method to reflect high and low flows. In this study, the high flow is reflected by annual 1-day maximum, 3-day maximum and 7-day maximum discharges calculated based on daily average discharges, and the low flow is reflected by annual 7-day minimum, 30-day minimum and 90-day minimum discharges calculated based on daily average discharges at three hydrological stations. These 6 indices we selected in this study are all contained in the indicators of hydrologic alterations (IHA) [81] to reflect the changes in flow regimes and have been widely used in previous studies [33,82,83].
To grasp the changes in future floods more accurately, the copula-based method is adopted to analyze future flood frequency. Before frequency analysis can be performed, flood events need to be identified. In this study, the peak over threshold (POT) method [84] is adopted. POT extracts flood events by setting a threshold value that leads to a specific number of flood events per year. With obvious advantages compared with the annual maximum flood series (AMF), POT can identify more flood events in wet years, ignore minor and insignificant floods in dry years, and provide more flood characteristics (such as flood volume, duration and lag time between adjacent events) [32,85]. This study adopted POT3 with three events per year, and the details of POT can be found in the studies by Lang et al. [84] and Mediero et al. [85].
Compared with classical bivariate frequency, the copula-based method performs marginal distribution and multivariate dependence modeling separately, providing greater flexibility in selecting marginal and joint probability functions [86,87]. In the field of bivariate frequency analysis, the copula function has been the dominant method and is widely used in flood frequency analysis. The copula function uses the following function to link the univariable marginal distribution functions to the multivariate probability distribution [88]: where F 1 (x 1 ), F 2 (x 2 ), . . . , F n (x n ) are marginal distributions of the random vector (x 1 , x 2 , . . . , x n ). A detailed description of copula functions can be found in the studies by Nelsen [89] and Salvadori et al. [90]. In this study, the flood peak and volume are selected as the marginal variables, and the candidate marginal distribution functions include the generalized extreme value (GEV), inverse Gaussian (IG), Birnbaum Saunders (BS), lognormal (LN), gamma (G), loglogistic (LL), Nakagami (Na), Rician (Ri), normal (N), logistic (L), t-location scale (T), Weibull (W), extreme value (EV), Rayleigh (Ra), exponential (Exp) and generalized Pareto (GP) distributions. The candidate copula functions include elliptical copulas (Gaussian), Archimedean copulas (Clatyon, Gumbel and Joe), extreme value copulas (Galambos), and mixture copulas (Plackett and Farlie-Gumbel-Morgenstren). All parameters were estimated using a hybrid-evolution Markov chain Monte Carlo (MCMC) approach within a Bayesian framework, and flood frequency analysis was performed by the Multivariate Copula Analysis Toolbox (MvCAT) [91,92]. A detailed mathematical description of marginal distribution functions, copula functions and joint return periods can be found in the study by Sadegh et al. [91].

Statistical Downscaling Results of Climate Variables
Using a combination of quantile mapping and machine learning, we performed spatial downscaling and bias correction for the coarse outputs of GCMs, and the results are shown in Figure 3. The comparison (Figure 3a) between the areal mean value of precipitation projections and observations showed that the downscaling method can satisfactorily correct the bias, and the future precipitation shows an increasing trend under the four SSP-RCPs. The fossil energy-developed SSP585 had the highest precipitation increase rate of 5.8 mm per year (0.33%/year), and the sustainably developed SSP126 had the second highest increase rate of 4.9 mm per year (0.28%/year) followed by the middle-of-the-road-developed SSP245 with an increase rate of 4.1 mm per year (0.23%/year), while the regional-rivalry-developed SSP370 had the lowest rate of increase of 2.3 mm per year (0.13%/year). In the last decade of the 21st century, annual precipitation is projected to increase by 15 (Figure 3b) shows that there are two peaks in a year according to the historical precipitation, which occur in June and August. In the future, two peaks still exist, with the second peak projected to increase greatly and almost surpass the first peak, which potentially indicates greater flood risks in the future. The projected annual mean temperature is shown in Figure 3c, and the results indicate that the future temperature shows different degrees of rise under the four SSP-RCPs. The sustainable development SSP126 has the slowest rise rate of 0.12 • C per decade, and continued warming will be effectively contained until 2070, subsequently showing a slight downward trend. In the last decade of the 21st century, the temperature under SSP126 will eventually increase by 0.96°C. The temperature under the other three SSP-RCPs will continue to rise steadily, and the rise rate will accelerate as greenhouse gas emissions increase with a rise rate of 0.16°C per decade under SSP245, a rise rate of 0.19°C per decade under SSP370 and 0.24°C per decade under SSP585. In the last decade of the 21st century, the annual mean temperature will increase 1.42°C under SSP245, 1.83°C under SSP370 and 2.34°C under SSP585 in the study area. In the case of constant land use, the gradual increase in temperature and precipitation in the future may increase evapotranspiration to some extent and eventually affect runoff. Since the forest coverage rate in the Jiulong River basin is nearly 57% and the evapotranspiration is largely contributed by vegetation, it is necessary to consider future land use or vegetation changes.

Projected Future Land Use Changes and Leaf Area Index under Various SSP-RCPs
Using the CA-Markov model with land use changes from LUH2 as constraints, this study projected decadal land use maps with a spatial resolution of 100 m × 100 m from 2020 to 2090. We summarize the areal mean land use changes of forest, cropland and urban areas in Figure 4. The historical land use map in 2010 and the projected land use maps in 2090 under four scenarios are shown in Figure 5, and other projected land use maps are provided in the Supplement (see Figures S1-S4). The results show that strict land use regulations and afforestation under SSP126 will lead to a continued increase in forest growth, with the ratio of forest coverage increasing by 10.4% in 2090 compared with that in 2010 and a decreasing trend of croplands, with the ratio decreasing by 10.7% compared with that in 2010. Under the other three SSP-RCPs, the ratio of forest will initially decrease and then increase, while the ratio of croplands will initially increase and then slowly decrease. With the lack of land use regulation, SSP370 shows the greatest decrease in forest and increase in croplands in 2050. Regarding urbanization, the SSP370 scenario will essentially

Projected Future Land Use Changes and Leaf Area Index under Various SSP-RCPs
Using the CA-Markov model with land use changes from LUH2 as constraints, this study projected decadal land use maps with a spatial resolution of 100 m × 100 m from 2020 to 2090. We summarize the areal mean land use changes of forest, cropland and urban areas in Figure 4. The historical land use map in 2010 and the projected land use maps in 2090 under four scenarios are shown in Figure 5, and other projected land use maps are provided in the Supplement (see Figures S1-S4). The results show that strict land use regulations and afforestation under SSP126 will lead to a continued increase in forest growth, with the ratio of forest coverage increasing by 10.4% in 2090 compared with that in 2010 and a decreasing trend of croplands, with the ratio decreasing by 10.7% compared with that in 2010. Under the other three SSP-RCPs, the ratio of forest will initially decrease and then increase, while the ratio of croplands will initially increase and then slowly decrease. With the lack of land use regulation, SSP370 shows the greatest decrease in forest and increase in croplands in 2050. Regarding urbanization, the SSP370 scenario will essentially maintain the current urbanization rate, while other SSP-RCPs will show a declining trend after 2050.  With projected land use maps and climate forcing data, the Biome-BGC model was used in this study to simulate future leaf area index (LAI) changes. First, the half-monthly remote sensing-based LAI was used to calibrate the Biome-BGC model from January 1985  With projected land use maps and climate forcing data, the Biome-BGC model was used in this study to simulate future leaf area index (LAI) changes. First, the half-monthly remote sensing-based LAI was used to calibrate the Biome-BGC model from January 1985 to December 1989 for each type of land use and to validate the model from January 1990 With projected land use maps and climate forcing data, the Biome-BGC model was used in this study to simulate future leaf area index (LAI) changes. First, the half-monthly remote sensing-based LAI was used to calibrate the Biome-BGC model from January 1985 to December 1989 for each type of land use and to validate the model from January 1990 to December 1994. The results ( Figure 6) mainly demonstrate the model performance for the LAI simulation of forest (including dense forest and sparse forest), cropland (including paddy field and upland field) and other vegetation types (including grasses and shrubs). The calibration and validation results showed that the model can basically capture the periodic changes in LAI with correlation coefficients greater than 0.4 for each type of land use, although certain discrepancies remain between the simulation and the remote-sensing LAI, which may be caused by the lack of detailed data of some local vegetation species and we must ignore these local species. With the calibrated Biome-BGC model, the future daily LAI under the four SSP-RCPs was projected.
Atmosphere 2022, 13, x FOR PEER REVIEW 12 of 22 the LAI simulation of forest (including dense forest and sparse forest), cropland (including paddy field and upland field) and other vegetation types (including grasses and shrubs). The calibration and validation results showed that the model can basically capture the periodic changes in LAI with correlation coefficients greater than 0.4 for each type of land use, although certain discrepancies remain between the simulation and the remote-sensing LAI, which may be caused by the lack of detailed data of some local vegetation species and we must ignore these local species. With the calibrated Biome-BGC model, the future daily LAI under the four SSP-RCPs was projected.

Future Changes in High and Low Flows
In this study, the GBHM was used to project future hydrological processes under the four SSP-RCPs with future climate forcing, land use maps and LAI projections, and the calibration and validation of the GBHM in the Jiulong River basin can be found in our previous study by Lu et al. [33]. The previous study has proven that the GBHM can accurately reflect the natural hydrological process in the study area and the calibration and validation results are shown in Table 2. PBIAS is the percent bias between the simulation and observation and can be used to evaluate the overall model performance [93]. NSE refers to the Nash and Sutcliffe efficiency coefficient and is widely used to evaluate the performance of a hydrological model [94]. LNSE is the logarithmic-transformed Nash and Sutcliffe efficiency coefficient and is used to evaluate the capability of low flow simulation of a hydrological model [95]. A hydrological model can be regarded as a satisfactory model when PBIAS is between −10% and 10% and NSE and LNSE are greater than 0.7

Future Changes in High and Low Flows
In this study, the GBHM was used to project future hydrological processes under the four SSP-RCPs with future climate forcing, land use maps and LAI projections, and the calibration and validation of the GBHM in the Jiulong River basin can be found in our previous study by Lu et al. [33]. The previous study has proven that the GBHM can accurately reflect the natural hydrological process in the study area and the calibration and validation results are shown in Table 2. PBIAS is the percent bias between the simulation and observation and can be used to evaluate the overall model performance [93]. NSE refers to the Nash and Sutcliffe efficiency coefficient and is widely used to evaluate the performance of a hydrological model [94]. LNSE is the logarithmic-transformed Nash and Sutcliffe efficiency coefficient and is used to evaluate the capability of low flow simulation of a hydrological model [95]. A hydrological model can be regarded as a satisfactory model when PBIAS is between −10% and 10% and NSE and LNSE are greater than 0.7 according to the study by Moriasi et al. [96]. The detailed calculation equations can be found in the study by Lu et al. [33].

Future Changes in High Flows
To make a preliminary assessment of future flood situations, we calculated the annual 1-day maximum, 3-day maximum and 7-day maximum discharges at three hydrological stations, with the results at the three stations shown in Figure 7. The annual 1-day maximum discharge showed the most variability, with the three stations increasing by an average of 133% under SSP126, 73% under SSP245, 68% under SSP370 and 67% under SSP585 in the last decade of the 21st century compared with that in the historical period. Annual 3-day maximum discharges had relatively less variability with the three stations increasing by an average of 141% under SSP126, 70% under SSP245, 72% under SSP370 and 64% under SSP585 in the last decade of the 21st century. The annual 7-day maximum discharges had the least variability, with the three stations increasing by an average of 144% under SSP126, 69% under SSP245, 87% under SSP370 and 64% under SSP585 in the last decade of the 21st century. The results demonstrate that the annual 1-day maximum, 3-day maximum and 7-day maximum discharges all reveal increasing trends in the future under all scenarios, which indicates that there may be more frequent floods with larger magnitudes in the future.
FOR PEER REVIEW 14 of 22 and 64% under SSP585 in the last decade of the 21st century. The annual 7-day maximum discharges had the least variability, with the three stations increasing by an average of 144% under SSP126, 69% under SSP245, 87% under SSP370 and 64% under SSP585 in the last decade of the 21st century. The results demonstrate that the annual 1-day maximum, 3-day maximum and 7-day maximum discharges all reveal increasing trends in the future under all scenarios, which indicates that there may be more frequent floods with larger magnitudes in the future. With daily projected runoff, we used the POT3 method to extract the flood events and then used the copula-based method to simulate the bivariate joint distribution of peak-volume combinations for the historical 30-year flood, 50-year flood and 100-year flood at the three stations under the four SSP-RCPs. To display the change in flood frequency more directly, the return period was adopted. A smaller return period represents the greater likelihood a given flood is to occur with a greater exceedance frequency, and vice versa. The results shown in Table 4 indicate that the frequency of floods of different magnitudes are projected to increase under the four future scenarios. The comparison among the four SSP-RCPs shows that the future described in SSP370 and SSP585 is projected to face more frequent floods, while the situation under SSP126 is projected to encounter the lowest increase in flood frequency, with flood frequency projections under SSP245 being in the middle. At Zhengdian station, the frequency of the historical 100-year flood is projected to increase by 128% under SSP126, 128% under SSP245, 184% under SSP370, and 178% under SSP585 in the last three decades of the 21st century. At Zhangping station, the frequency of the historical 100-year flood is projected to increase by 86% under SSP126, 120% under SSP245, 239% under SSP370, and 263% under SSP585. At Punan station, the frequency of the historical 100-year flood is projected to increase by 151% under SSP126, 165% under SSP245, 215% under SSP370, and 158% under SSP585. Increasingly frequent floods in the future may cause serious damage to the lives and properties of local people in the study area, and effective mitigation and adaptation measures should be proposed to face future increasing flood risks. With daily projected runoff, we used the POT3 method to extract the flood events and then used the copula-based method to simulate the bivariate joint distribution of peakvolume combinations for the historical 30-year flood, 50-year flood and 100-year flood at the three stations under the four SSP-RCPs. To display the change in flood frequency more directly, the return period was adopted. A smaller return period represents the greater likelihood a given flood is to occur with a greater exceedance frequency, and vice versa. The results shown in Table 4 indicate that the frequency of floods of different magnitudes are projected to increase under the four future scenarios. The comparison among the four SSP-RCPs shows that the future described in SSP370 and SSP585 is projected to face more frequent floods, while the situation under SSP126 is projected to encounter the lowest increase in flood frequency, with flood frequency projections under SSP245 being in the middle. At Zhengdian station, the frequency of the historical 100-year flood is projected to increase by 128% under SSP126, 128% under SSP245, 184% under SSP370, and 178% under SSP585 in the last three decades of the 21st century. At Zhangping station, the frequency of the historical 100-year flood is projected to increase by 86% under SSP126, 120% under SSP245, 239% under SSP370, and 263% under SSP585. At Punan station, the frequency of the historical 100-year flood is projected to increase by 151% under SSP126, 165% under SSP245, 215% under SSP370, and 158% under SSP585. Increasingly frequent floods in the future may cause serious damage to the lives and properties of local people in the study area, and effective mitigation and adaptation measures should be proposed to face future increasing flood risks.

Future Changes in Low Flows
To make a preliminary assessment of future hydrological drought situations, we calculated the annual 7-day minimum, 30-day minimum and 90-day minimum discharges at three hydrological stations, and the results are shown in Figure 8. For changes in low flows, the results show that there are slightly increasing trends in the 7-day minimum, 30-day minimum and 90-day minimum discharges in the future at the three stations under the four SSP-RCPs compared with those in the historical periods. The annual 7-day minimum discharge had the least variability, with the three stations increasing by an average of 39% under SSP126, 19% under SSP245, 25% under SSP370 and 0.26% under SSP585 in the last decade of the 21st century compared to that in the historical period. Annual 30-day minimum discharges have relatively greater variability with the three stations increasing by an average of 41% under SSP126, 18% under SSP245, 24% under SSP370, and 31% under SSP585 in the last decade of the 21st century. The annual 90-day minimum discharges had the most variability, with the three stations increasing by an average of 36% under SSP126, 13% under SSP245, 25% under SSP370, and 29% under SSP585 in the last decade of the 21st century. This indicates that future hydrological droughts may be alleviated to some extent, although droughts have not been a severe problem in the Jiulong River basin.
x FOR PEER REVIEW 16 of 22

Impact of Climate and Land Use Changes on Future High and Low Flows
The above results have demonstrated a wetter and warmer future in this coastal basin and analyzed potential changes in future runoff and flood frequency under the combined impacts of climate change and land use change. However, the respective impacts of climate change and land use change need to be further discussed. Here, we replaced the historical 30-year climate forcing data with the projected 30-year (2071-2100) future climate data under SSP370 and retained the historical land use (referred to as the newly designed scenario). The LAI was simulated by the calibrated Biome-BGC model, and the hydrological process was projected using the calibrated GBHM model. Since SSP370 and the newly designed scenario had the same precipitation in 2071-2100 with the only difference being in land use change, the comparison between these two scenarios could show the impact of land use change. Similarly, the only difference observed between the historical situation and the newly designed scenario was the change of climate forcing data, and the comparison between these two scenarios could quantify the impact of climate change. The annual 1-day maximum discharge was selected to represent high flow, and the annual 7-day minimum discharge was chosen to represent low flow.
For high flow, the comparison results showed that the mean value of annual 1-day maximum discharge under the newly designed scenario is projected to increase by 106% at Zhengdian station, 100% at Zhangping station and 48% at Punan station compared to that in the historical period , while that under SSP370 is projected to increase by 128% at Zhengdian station, 103% at Zhangping station and 51% at Punan station. It can be concluded that climate change is a dominating factor for the future increase in high flows with larger precipitation during the wet season. The comparison between the newly designed scenario and SSP370 shows that continuous deforestation and cropland exten-

Impact of Climate and Land Use Changes on Future High and Low Flows
The above results have demonstrated a wetter and warmer future in this coastal basin and analyzed potential changes in future runoff and flood frequency under the combined impacts of climate change and land use change. However, the respective impacts of climate change and land use change need to be further discussed. Here, we replaced the historical 30-year climate forcing data with the projected 30-year (2071-2100) future climate data under SSP370 and retained the historical land use (referred to as the newly designed scenario). The LAI was simulated by the calibrated Biome-BGC model, and the hydrological process was projected using the calibrated GBHM model. Since SSP370 and the newly designed scenario had the same precipitation in 2071-2100 with the only difference being in land use change, the comparison between these two scenarios could show the impact of land use change. Similarly, the only difference observed between the historical situation and the newly designed scenario was the change of climate forcing data, and the comparison between these two scenarios could quantify the impact of climate change. The annual 1-day maximum discharge was selected to represent high flow, and the annual 7-day minimum discharge was chosen to represent low flow.
For high flow, the comparison results showed that the mean value of annual 1-day maximum discharge under the newly designed scenario is projected to increase by 106% at Zhengdian station, 100% at Zhangping station and 48% at Punan station compared to that in the historical period , while that under SSP370 is projected to increase by 128% at Zhengdian station, 103% at Zhangping station and 51% at Punan station. It can be concluded that climate change is a dominating factor for the future increase in high flows with larger precipitation during the wet season. The comparison between the newly designed scenario and SSP370 shows that continuous deforestation and cropland extension under SSP370 will exacerbate this increase (from 106% to 128% at Zhengdian station, 100% to 103% at Zhangping station and 48% to 51% at Punan station) by affecting canopy interception, infiltration and evapotranspiration processes, which is consistent with previous studies [97][98][99]. For low flow, the results demonstrate that the mean value of annual 7-day minimum discharge under the newly designed scenario will increase by 8.5% at Zhengdian station, 8.4% at Zhangping station and 9.5% at Punan station, while that under SSP370 will increase by 27.2% at Zhengdian station, 31.3% at Zhangping station and 19.3% at Punan station. This indicates that future climate change in the Jiulong River basin may have less impact on low flows, while land use changes (deforestation and cropland extension) will cause an increase in low flows. The main reason for this phenomenon is that forest degradation reduces evapotranspiration and canopy interception, which leads to an increase in soil water content and groundwater recharge and ultimately causes an increase in low flow [1,100].

Impact of Dams on Futrue High and Low Flows
Since there are numerous dams existing in the Jiulong River basin, the possible impacts of dams on future high and low flows are desired points of discussion. However, since the GBHM model we used in this study did not have the reservoir operation module which requires further improvement, it was difficult to quantitatively analyze the current impact of dams thus we only quantitatively discussed the possible impacts here. The drainage area of Zhengdian station is mainly affected by annual flood-control dams according to the reservoir capacities, and the drainage area of Zhangping and Punan stations are mainly affected by small-scale hydropower dams (<50 MW). Our previous study [33] showed that the flood-control dams could decrease the high flows and increase the low flows at Zhengdian station by releasing water in dry seasons and restraining floods in wet seasons, while at Zhangping and Punan stations, the dominant hydropower dams could reduce low flows and had less obvious effects on high flows. This is due to these small dams being built solely to generate as much electricity as possible and may thus store more water up to a sufficient head before electricity generation can occur during the dry season. Under a future with increasing high flows and frequent floods, the flood-control dams may to some extent mitigate the growth of high flows and reduce flood risks in addition to increasing low flows at Zhengdian station. However, the numerous small-scale hydropower dams may have less effect on the reduction of high flows and mitigation of flood risks at Zhangping and Punan stations, and the slight increase in low flows due to climate and land use changes may be offset to some extent by the reduction effect from hydropower dams, thus low flows at these two stations may not change much in the future.

Suggestions for Future Water Management
According to the results found in this study, the Jiulong River basin faces a challenging future with more frequent floods, and a wetter and warmer climate will dominate the increase in high flows, which will be further exacerbated by forest degradation. Current water management systems will be severely challenged and need to be adjusted. Integrated water management is required in the future. First, as regional positive and sustainable land use changes can offset the impacts of climate change to some extent, we suggest that local authorities implement strict land management practices and vigorously promote afforestation. Second, we recommend that the flood protection standards of the infrastructure should be improved in the Jiulong River basin, which is in close proximity to several cities with large populations and thriving economies. Third, we suggest that the management department should improve the operation rules of numerous existing reservoirs to adapt to more severe flood situations in the future.

Conclusions
With a combined downscaling quantile mapping-machine learning method, the future climate forcing data from GCMs in CMIP6 were bias corrected and spatially downscaled. Using the CA-Markov model with future land use changes from LUH2 as constraints, decadal future land use maps with high spatial resolution were projected. Based on the climate forcing and land use projections, the Biome-BGC model was applied to project future LAI data, and subsequently, a calibrated hydrological model (GBHM) was used to simulate the future hydrological process under four SSP-RCPs during the period from 2015 to 2100. The main conclusions were as follows: (1) Under global warming, precipitation and air temperature will continue to rise under all SSP-RCPs. In the last decade of the 21st century, precipitation is projected to increase by 15.5% under SSP126, 22.2% under SSP245, 14.6% under SSP370 and 24.9% under SSP585, and the annual mean temperature will increase by 0.96°C under SSP126, 1.42°C under SSP245, 1.83°C under SSP370 and 2.34°C under SSP585 in the study area.
(2) In the future, land use changes will have different development pathways determined by the initial designs in various SSP-RCPs. Under greenest roads (SSP126), the forest area will continue to rise with the ratio increasing by 10.4% in 2090, and the cropland decreasing by 10.7% compared with that in 2010. Under the other three SSP-RCPs, the ratio of forest will initially decrease and then increase, while the ratio of croplands will initially increase and then slowly decrease.
(3) Under the combined impacts of climate and land use change, the total available water resources will increase due to increasing precipitation, and the high flow and low flow will both increase at three stations under the four SSP-RCPs. In addition, the Jiulong River basin will face more frequent floods in the future. In the last three decades of the 21st century, the station-average frequency of a historical 100-year flood will increase by 122% under the most optimistic of scenarios (SSP126), and increase by 213% under the scenario of greatest regional rivalry (SSP370).
(4) Future warmer and wetter climates will lead to a great increase in high flows, and deforestation will exacerbate this increase to some extent, while the increase in low flows mainly caused by land use changes (deforestation and cropland extension), and the impact of climate change will be relatively smaller.
(5) To tackle the challenges of future water management presented by climate change, we recommend an integrated water management system, including afforestation, raised protection standards of infrastructure and improvement of existing operation rules of numerous reservoirs.  Data Availability Statement: GIMMS LAI3g is available at https://ecocast.arc.nasa.gov/data/pub/ gimms/3g.v1/, accessed on 15 January 2022. The GCM data can be directly downloaded at https: //esgf-node.llnl.gov/search/cmip6/, accessed on 15 January 2022. The LUH2 dataset is available at https://luh.umd.edu/data.shtml, accessed on 15 January 2022. The other datasets analyzed during the current study are available from the corresponding author upon reasonable request.