The NPP-Based Composite Indicator for Assessing the Variations of Water Provision Services at the National Scale

Water provision (WP) is an important service of the terrestrial ecosystem, which contributes to water availability for consumptive use and in situ water supply, sustains the production or flows of multiple ecosystem services (ES). Spatially explicit mapping of WP is critical for incorporating the ES concept into the decision-making processes of land-use and ecological conservation planning. Traditionally, regional complexes hydrological process models were simplified and used for mapping WP of the ecosystem at broad scales, but this approach is significantly limited by data accessibility and difficulty validating the results. To fill the gap, an NPP-based composite indicator model that simulates WP by multiplying NPP and its variations with the soil infiltration capacity factor, annual precipitation and the slope of the land surface is proposed in this paper. These parameters are chosen to map WP because they are closely related to hydrological processes. The model results were validated using observed runoff data of the eleven river basins in China. We then applied this approach to analyze the spatiotemporal changes of WP in China from 2000 to 2013. The results show that: (1) the average value of WP was lowest in the Northwest Arid Area ecoregions while the highest value of WP was in the South China ecoregion. (2) The linear trend of WP in the Loess Plateau and Hengduan Mountains ecoregions were increased while decreased in the other nine ecoregions. (3) The WP in the north of the Qinghai-Tibet Plateau presented a significant decrease trend mostly because the land cover change (e.g., grassland convert into dessert) and decreasing precipitation; the decreasing of the WP in Yunan-Guizhou Plateau are because the farmland convert into settlement land and the significant decrease of precipitation and significantly increase of temperature; the significant increase of the WP in Northeast China are because the increasing of forest and farmland, the grassland and wetland convert into farmland and forest, and the significant decrease of temperature and increase of precipitation; Although the increase of precipitation has played an important role in promoting WP, the significant increase of WP in the Loess Plateau was mainly due to the farmland convert into forest and grassland ecosystem types. The indicator explored by this research is benefiting for revealing the variations of WP under different land-use change and climate change, and informed the decision-making process of land-use policy or conservation planning at data-scarce regions or broaden spatial scales.

to 2013. Finally, we discuss the applicability of the model to assess WP across spatial scales, and the potential WP conservation policy implication of the findings.

Study Area
China, with the land area of 9.6 million square kilometers, is character by diverse geomorphic and climate types, populations (1.395 billion by 2018) and social-economic development levels. The main ecosystem types of China are include grassland (31.37%), forest (24.06%), farmland (19.6%), desert (13.41%), other (6.49%), wetland (3.54%) and settlement (1.54%) (Figure 1). Along latitude gradient from the south to the north, its tropical zone, warm temperate zone and the cold temperate zone [27,28]. China has rich river network systems. The well-known seven major river systems in China include the Pearl River system, Yangtze River system, Yellow River system, Huaihe River system, Liaohe River system, Haihe River system and Songhuajiang River system. These seven major river basins cover around 64 percent of the country's total land area. Based on precipitation, temperature, altitude, vegetation index and geomorphological conditions, the land territory of China were divided into eleven ecoregions [29]: the Northeast China (1), Northwest arid area (2), Inner Mongolia plateau (3), Qinghai-Tibet Plateau (4), the Loess Plateau (5), the northern China (6), the Hengduan Mountains (7), the Sichuan and Chongqing region (8), the Middle and lower reaches of the Yangtze River (9), the Yunan-Guizhou Plateau (10) and the South China (11) (Figure 1).
Water 2019, 11, x FOR PEER REVIEW 3 of 13

Study Area
China, with the land area of 9.6 million square kilometers, is character by diverse geomorphic and climate types, populations (1.395 billion by 2018) and social-economic development levels. The main ecosystem types of China are include grassland (31.37%), forest (24.06%), farmland (19.6%), desert (13.41%), other (6.49%), wetland (3.54%) and settlement (1.54%) (Figure 1). Along latitude gradient from the south to the north, its tropical zone, warm temperate zone and the cold temperate zone [27,28]. China has rich river network systems. The well-known seven major river systems in China include the Pearl River system, Yangtze River system, Yellow River system, Huaihe River system, Liaohe River system, Haihe River system and Songhuajiang River system. These seven major river basins cover around 64 percent of the country's total land area. Based on precipitation, temperature, altitude, vegetation index and geomorphological conditions, the land territory of China were divided into eleven ecoregions [29]: the Northeast China (1), Northwest arid area (2), Inner Mongolia plateau (3), Qinghai-Tibet Plateau (4), the Loess Plateau (5), the northern China (6), the Hengduan Mountains (7), the Sichuan and Chongqing region (8), the Middle and lower reaches of the Yangtze River (9), the Yunan-Guizhou Plateau (10) and the South China (11) (Figure 1).

The NPP-Based Model for Mapping WP
The net primary production (NPP) of vegetation combined with a variety of other ecosystem attributes are used as a composite indicator model for mapping the WP service at a large scale [23,26] (Equation (1)): where WP refers to the capability of the water provision service of the ecosystem. The NPP is the annual net primary production of vegetation. The VCnpp describes the standard deviation of the NPP

The NPP-Based Model for Mapping WP
The net primary production (NPP) of vegetation combined with a variety of other ecosystem attributes are used as a composite indicator model for mapping the WP service at a large scale [23,26] (Equation (1)): where WP refers to the capability of the water provision service of the ecosystem. The NPP is the annual net primary production of vegetation. The VC npp describes the standard deviation of the NPP within a year. We used the terrestrial Carnegie Ames-Stanford Approach (CASA) to estimate the NPP of ecosystems [30] (Equation (2)). The VC npp is the standard deviation of the monthly NPP within a year. Soil infiltration capacity factor F sic is derived from the USDA (United States Department of Agriculture) soil texture dataset, with 13 categories of soil particles [31]. In ascending order, they are clay (heavy), silty clay, clay, silty clay loam, clay loam, silt, silt loam, sandy clay, loam, sandy clay loam, sandy loam, loamy sand, and sand. Each category is assigned a value of 1/13, 2/13, . . . , and 13/13, following the same order. Annual precipitation map F pre was interpolated from 837 meteorological station data in China using the ordinary Kriging method in the Geostatistical Analyst modular of the ArcGIS 10.2 software (ESRI) (Redlands, CA, USA). The F slo refers to the slope of the land surface. The terrestrial Carnegie Ames-Stanford Approach (CASA) was used to estimate the NPP of vegetation. The CASA model advocates determining the NPP of vegetation as the product of the modulated absorbed photosynthetically active radiation (APAR) and the light use efficiency (ε) factor [30,32]: where NPP(x,t) describes the net primary production of location x at time t, APAR refers to the canopy-absorbed incident solar radiation (MJ × m −2 ), and the parameter ε(x,t) is the light use efficiency (g × C·× MJ −1 ). Except for NPP, all the parameters of VC npp , F slo , and F pre factors are standardized between 0 and 1 by the maximum and minimum values of their values from the year 2000 to 2013. The units of the WP are dimensionless. In this study, the value of WP in China is ranged from 0 to 466.51. The zero WP regions are located at the desert areas of the Northwest arid ecoregion, where no vegetation cover is distributed and the annual precipitation is zero. The largest WP regions are located in the tropical monsoon rainforest of Hainan Island of South China ecoregion where the NPP and annual precipitation is the largest in China. All the spatial analysis (e.g., interpolation and raster calculation) was performed in ArcGIS 10.2 software (ESRI).

Other Analysis Methods
The linear regression trend approach has been widely used in trend detection for vegetation change and climate variations in a continuous period [33]. In this paper, we used linear regression to analyze the changing trend of WP, precipitation, and temperature from 2000 to 2013. The classification of linear trend is determined by the value of the linear trend and its significance level. It is classified into five categories: significantly decrease (P < 0.05 and linear trends a < 0), no changes (a = 0), no significant decrease (P > 0.05 and a < 0), no significant increase (P > 0.05 and a > 0), and significantly increase (P < 0.05 and a > 0). Pearson's correlation analyses were performed using the function of [r,p] = corr (X,Y) in MATLAB Version R2010b (Natick, MA, USA). [33], which can return both the correlation coefficients r and p-value (significant test) simultaneously, where r and p represent the correlation coefficient matrix and the matrix of significant test results and X and Y are the respective vectors of WP and climate variation (e.g., annual precipitation or temperature) characteristics from 2000 to 2013. The significance of correlations between climate factors and annual WP are classified into five categories: significantly negative (P < 0.05 and the Person correlation coefficient (R) < 0), no significant negative (P > 0.05 and R < 0), no significant positive (P > 0.05 and R > 0), significant positive (P < 0.05 and R > 0), and no relations (R = 0).

Data Sources
The data used for mapping WP were as follows: daily meteorological data (solar radiation, precipitation, and temperature) from 2000 to 2013 were retrieved from the China Meteorological Data Sharing Service System (https://data.cma.cn/site/index.html) (Figure 1). Data at 837 stations within China were used to produce meteorological raster maps at 250m resolution using the ordinary Kriging method in the Geostatistical Analyst modular of the ArcGIS 10.2 (ESRI) software. The 250 m MODIS NDVI data acquired from the Atmosphere Archive and Distribution System were used to calculate the annual total NPP (https://ladsweb.modaps.eosdis.nasa.gov/). Spatial distribution data of terrestrial ecosystem types in China in 2000 and 2010 were collected from the Resource and Environment Data Cloud Platform (http://www.resdc.cn/). Topographical parameters were derived from the Shuttle Radar Topography Mission (STRM) digital elevation data at 90 m resolution (http://srtm.csi.cgiar.org/). The soil properties data was provided by the Cold and Arid Regions Sciences Data Center in Lanzhou (http://westdc.westgis.ac.cn). All the raster maps data were interpolated or resampled to 250 m resolution in the ArcGIS10.2 software (ESRI) before being input into the WP model for further analysis.

Evaluate of the Model Result
Generally, the validation of the processes-based model at the local scale is relatively simple to use the observation data. However, the validation of the regional model is very difficult because they usually applications based on parameterization rules developed at the local scale are performed at the regional scale where only gauge data are available for comparison [12]. In addition, due to data and parameter uncertainty and the fact that the distributed model cannot be calibrated at the regional scale, the model results remain uncertain [12]. Therefore, the accurate assessment of "available water" is difficult, especially in data-scarce regions [5]. Although many scholars believe that it is impossible to verify numerical models applied to large scales because the complexity and heterogeneity of geographical and ecological environment bring many uncertainties for verification and validation [10,23,34]. However, in this study, the runoff of the eleven highly representative sub-basins of the seven major river systems in China was used to evaluate the WP model results from 2000 to 2013 ( Figure 2). The Pearson correlation coefficient between the WP and the runoff of the sub-basin were all higher than the significance test of 0.05, which indicates that the explored WP model has good accuracy and availability on a large scale. Therefore, we believe that the NPP-based model is a suitable approach for mapping WP and can be used in the further analysis of spatial patterns and linear trend assessment processes ( Figure 2).

Evaluate of the Model Result
Generally, the validation of the processes-based model at the local scale is relatively simple to use the observation data. However, the validation of the regional model is very difficult because they usually applications based on parameterization rules developed at the local scale are performed at the regional scale where only gauge data are available for comparison [12]. In addition, due to data and parameter uncertainty and the fact that the distributed model cannot be calibrated at the regional scale, the model results remain uncertain [12]. Therefore, the accurate assessment of "available water" is difficult, especially in data-scarce regions [5]. Although many scholars believe that it is impossible to verify numerical models applied to large scales because the complexity and heterogeneity of geographical and ecological environment bring many uncertainties for verification and validation [10,23,34]. However, in this study, the runoff of the eleven highly representative sub-basins of the seven major river systems in China was used to evaluate the WP model results from 2000 to 2013 ( Figure 2). The Pearson correlation coefficient between the WP and the runoff of the sub-basin were all higher than the significance test of 0.05, which indicates that the explored WP model has good accuracy and availability on a large scale. Therefore, we believe that the NPP-based model is a suitable approach for mapping WP and can be used in the further analysis of spatial patterns and linear trend assessment processes ( Figure 2).

Figure 2.
Spatial distributions of the eleven sub-basins for validation. R represents the correlation coefficients between WP and runoff in each sub-basin. The double asterisk (**) indicates that the significance level of correlation coefficient (R) is less than 0.01 and the asterisk (*) represents the significance level of R is greater than 0.01 and less than 0.05. Figure 2. Spatial distributions of the eleven sub-basins for validation. R represents the correlation coefficients between WP and runoff in each sub-basin. The double asterisk (**) indicates that the significance level of correlation coefficient (R) is less than 0.01 and the asterisk (*) represents the significance level of R is greater than 0.01 and less than 0.05.   (Figure 3b). Although the annual WP in China from 2000 to 2013 slightly declined (−1.08 unit per year), the areas with increased WP (62.51%) were larger than the areas with decreased WP (36.49%). Furthermore, the significantly increased regions (11.43%) were larger than the significantly decreased regions (3.79%) ( Table 1).  (Figure 3b). Although the annual WP in China from 2000 to 2013 slightly declined (−1.08 unit per year), the areas with increased WP (62.51%) were larger than the areas with decreased WP (36.49%). Furthermore, the significantly increased regions (11.43%) were larger than the significantly decreased regions (3.79%) ( Table 1). The significantly changed regions also had distinct spatial heterogeneity (Figure 3b). Among the significantly increased region, the WP in Loess Plateau ecoregions increased most (39.27%); the Northeast ecoregions come second, 29.27%; the other nine ecoregions with slightly increase. Among the significantly decreased regions, the WP in Yunnan-Guizhou Plateau ecoregion decreased the most (25.01% regions decreased) ( Table 1). The significantly increased regions centered in the west part of the Northeast and the northern of Loess Plateau ecoregions, the significantly decreased regions gathered in the center of the Yunnan-Guizhou Plateau ecoregion, while the other significantly changed regions scattered across the China mainland.

Key Influence Factors of the Variations in WP
The amount of WP is limited by the precipitation in time and space and the spatial patterns of WP ecosystem are highly related to the spatial distributions of precipitation, vegetation, landforms, human-induced changes in land use and consuming activities [23,[35][36][37]. The annual variations of WP service were primarily influenced by land use changing and climatic variations [23,35]. To investigate the quantitative relationship between WP and climate change (precipitation and temperature), we performed a spatial correlation between the time-series maps ( Figure 4). The spatial heterogeneity of the linear trend of temperature is very high. The significantly decreased region of temperature is more located at the Northwest arid area while the significantly increased area is more located at the Qinghai-Tibet Plateau and the Yunan-Guizhou Plateau. The linear correlation between annual temperature and WP is significantly negative in Northeast China and the Yunan-Guizhou Plateau. The significantly increased regions of precipitation are more located in Northeast China while the significantly decrease regions is more located at the Qinghai-Tibet Plateau and the Yunan-Guizhou Plateau. Nearly all of the area of precipitation has a positive correlation with WP ( Figure 4).
Besides, the transformation of ecosystems from 2000 to 2010 was also used to reveal the effects of land use on WP (Table 2). From the whole China scale, the area of farmland and grassland increased while the critical transformation characteristics of the ecosystem are the convert of farmland into the settlement ecosystem. In the significantly increase region of WP, the area of farmland and forest is increased that is because of the convert of grassland, wetland, and desert into farmland and forest ecosystem. In the Qinghai-Tibet Plateau, the typical feature of ecosystem variation is the grassland convert into dessert and settlement while the dessert converts into the wetland. In the Loess Plateau, the main transformation characteristics are the convert of farmland into forest and grassland ecosystems. In Northeast China, the increasing of forest and farmland are caused by the convert of grassland and wetland ecosystem. In the in Yunan-Guizhou Plateau, the decreasing farmland is because of the convert of farmland into the settlement ecosystem. The coastal significantly decreasing areas of the Middle and lower reaches of the Yangtze River and South China are located in the giant city of Shanghai and Guangzhou. In the process of rapid urbanization in these two cities, farmland and forestland are gradually eroded by urban land, resulting in a significant decrease of WP (Figures 1 and 3b). Besides, the transformation of ecosystems from 2000 to 2010 was also used to reveal the effects of land use on WP (Table 2). From the whole China scale, the area of farmland and grassland increased while the critical transformation characteristics of the ecosystem are the convert of farmland into the settlement ecosystem. In the significantly increase region of WP, the area of farmland and forest is increased that is because of the convert of grassland, wetland, and desert into farmland and forest ecosystem. In the Qinghai-Tibet Plateau, the typical feature of ecosystem variation is the grassland convert into dessert and settlement while the dessert converts into the wetland. In the Loess Plateau, the main transformation characteristics are the convert of farmland into forest and grassland ecosystems. In Northeast China, the increasing of forest and farmland are caused by the convert of grassland and wetland ecosystem. In the in Yunan-Guizhou Plateau, the decreasing farmland is because of the convert of farmland into the settlement ecosystem. The coastal significantly decreasing areas of the Middle and lower reaches of the Yangtze River and South China are located in the giant city of Shanghai and Guangzhou. In the process of rapid urbanization in these two cities, farmland and forestland are gradually eroded by urban land, resulting in a significant decrease of WP (Figures 1 and 3b).  In sum, the significant decrease of the WP in the north of the Qinghai-Tibet Plateau was due to the grassland convert into dessert and the dessert convert into wetland, and the significant decrease of precipitation; the decreasing of the WP in Yunan-Guizhou Plateau are because the farmland convert into settlement ecosystem, the significant decrease of precipitation and significantly increase in temperature; the significant increase of the WP in Northeast China are because the increasing of forest and farmland, the grassland and wetland convert into farmland and forest, the significantly decrease of temperature and increase of precipitation; Although the increase of precipitation has played an important role in promoting WP, the significant increase of WP in the Loess Plateau was mainly due to the farmland convert into forest and grassland ecosystem types. The significant decreasing of WP in coastal areas of the Middle and lower reaches of the Yangtze River and South China are because of the farmland and forestland gradually eroded by urban land during the urbanization process of Shanghai and Guangzhou (Figures 3b and 4, Table 2).

Limitations and Usefulness in Future Studies
Some limitations and usefulness can be concluded when using this method in further studies. Firstly, because the main indicator of this composite indicator method is NPP data, which contributes around 80% of the WP value [23]. Therefore, the spatial distribution of the hotspots for conserving WP will mostly be concentrated in forestland [23,25,26]. However, from the perspective of ecosystem diversity protection, the other types of ecosystem types also make important contributions to sustain WP services in given spatial scales. Therefore, in the specific designation of WP protection areas, diversity factors of different ecosystem types should be included in a conservation program at the specific ecological regions. For example, both incorporate the total amount of WP and the area proportions of each ecosystem types in the prioritization processes of WP conserving locations. Secondly, although the NPP-based composite indicator for mapping WP can reveal the relative rankings of WP in a specific region rather than the specific amounts of WP modeled by using process-based approaches, they are sufficient to meet the requirements of priority setting for conservation or land-use management [22][23][24][25][26]. Finally, as discussed in the above sections, many scholars believe that it is difficult to spatial explicitly accurately assess and verify the amount of water provide in each component of a hydrological process due to the complexity of the process, the scarcity of data, and the uncertain impact of numerous natural and human factors [10,12,23,34]. So the shortcoming of this easy to use NPP-based model was the validation of the results. Generally, the precipitation over land is divided into the green water and blue water resources [36,37]. The blue water is the runoff via groundwater and surface water, and the green water is the rainfall that intercepted by vegetation and eventually returns to the atmosphere as evapotranspiration or infiltrates into the soil [36,37]. The amount of blue water is usually less than green water and both are limited by the precipitation in time and space [37]. The patterns of WP in the landscape is highly affected by the spatial distributions of rainfall, vegetation, human-induced changes in land use and consuming activities [35][36][37][38]. The NPP-based composite indicator model of this article without dividend the WP into the green and blue water, which simulates WP by multiplying NPP and its variations with the soil infiltration capacity factor, annual precipitation and the slope of the land surface [23][24][25][26]. These parameters are chosen to map WP because they are closely related to hydrological processes and both integrate the green water-related parameter and blue water-related variables [23,36,37]. Through the coupling analysis of the spatial-temporal change pattern of WP services on a national scale and its driving factors with climate change and land-use change, it can be seen that the composite indicator method based on NPP can not only meet these requirements but also effectively reveal the driving factors of WP service changes. The results show that the NPP-based method for mapping WP is sensitive to runoff, precipitation, and land-use change. Although there is room for future improvement the NPP-based model to build composite indicators for assessing multiple ES, the opportunities of this method are the convenience of using it in multi-scales, such as local, regional and national policy decision making processes or in the data-scarce regions [10,26]. Furthermore, the sustainable use of WP for different purposes will not only need to consider the demand in water quantity, but also for water temperature and salinity, nutrient levels and other pollutants [39].The presented analysis only focuses on water quantity, but water quality is also an important factor affecting available freshwater resources to ecosystems and human use. Globally, there is a growing concern that water quality degradation, such as non-point source pollution from agricultural production [40], should be considered when evaluating WP capabilities [39]. It is possible to extend the NPP-based model to consider the impacts of different land-use types on nutrient loading and therefore WP capabilities. The environmental degradation resulting from socio-economic development is prevalent in many developing countries, this approach can be used for informing ecological conservation decisions in countries or regions facing similar problems.

Summary and Conclusions
The consistent and rapid growth of population and economic development in recent decades poses increasing environmental challenges, especially the water crisis [41,42]. To restore and sustain the fragile ecosystems, the Chinese government has proposed several ecological protection measures to curb the trend of environmental degradation [6,7,43]. These ecological protection measures have delineated key areas for water conservation services. However, due to the different starting times of these protected areas and the absence of strict legal protection, these protected areas have essentially faced or suffered serious environmental damage over time [6,43,44]. Therefore, both from the perspective of decision-making departments and scientific research, it is necessary to put forward practical and effective evaluation methods for the effect of rapid ecological protection areas suitable for large-scale applications and to formulate effective ecological compensation plans and ecological restoration and protection measures [45]. Among the existing evaluation methods, the ecological process-based method is the preferred long-term evaluation and detection method because it is easy to be accepted by experts and scholars. However, this method simplifies the application of model parameters established based on long-term data observation on a small scale to a large scale, and the results will inevitably bring some uncertainty at large scale, for example, the results usually lack spatial location information. Therefore, it is urgent to propose a simple, effective and relatively accurate evaluation methods to monitor changes in WP services in key ecoregions. The NPP-based composite indicator model of this article may meet these policy gaps.
In sum, the WP service is critical for economic development and environmental sustainability. In this study, an NPP-based composite indicator approach was proposed to map the WP on a broad scale, with China as a case study. The sensitive revelation of the model result of spatial patterns changing of WP service and its driving factors on a large scale indicates the practicability and reliability of the method proposed. The results show that the average value of WP was lowest in the Northwest Arid Area ecoregions while the highest value of WP was in the South China ecoregion. The linear trend of WP in the Loess Plateau and Hengduan Mountains ecoregions were increased while decreased in the other nine ecoregions. The significant variations of WP in the north of the Qinghai-Tibet Plateau, Yunan-Guizhou Plateau, Northeast China, and the Loess Plateau was mainly due to land use and climate-changing. The significant decreasing of WP in coastal areas of the Middle and lower reaches of the Yangtze River and South China are because of the urbanization process of Shanghai and Guangzhou. In brief, the indicator explored by this research is affecting and benefiting for revealing the variations of WP under different land use changing and climate scenario backgrounds, these findings can help inform national water conservation planning or swift assessing and longtime surveillance of WP variations in conservational areas.