Greening Trends of Southern China Confirmed by GRACE

As reported by the National Aeronautics and Space Administration (NASA), the world has been greening over the last two decades, with the highest greening occurring in China and India. The increasing vegetation will increase plant tissue accumulation and water storage capacity, and all of these variations will cause mass change. In this study, we found that the mass change related to greening in Southern China could be confirmed by Gravity Recovery and Climate Experiment (GRACE) observations. The mean mass change rate detected by GRACE is 6.7± 0.8 mm/yr in equivalent water height during 2003–2016 in our study region. This is consistent with the sum of vegetation tissue, soil water and groundwater change calculated using multi-source data. The vegetation accumulation is approximately 3.8 ± 1.3 mm/yr, which is the major contribution to region mass change. We also found that the change of water storage capacity related to vegetation can be detected by GRACE.


Introduction
Vegetation is essential for the subsistence of all life on Earth, as it provides oxygen through photosynthesis and consumes greenhouse gases to slow the rate of global warming. Vegetation change has a large influence on water resources, and this has been highly recognized by the Intergovernmental Panel on Climate Change (IPCC) [1]. Glacial ice is the largest source of freshwater on the Earth. However, it was reported in recent studies that ice sheets and high mountain glaciers are decreasing in area in response to global warming [2,3], which will affect the global trends of freshwater availability [4]. Surface runoff is a major component of the water cycle which affects land water storage change. Wei et al. [5] found that the variation of vegetation cover has a strong effect on the annual global average variation of runoff, and the greater the forest losses, the greater the annual runoff variation. According to Bhanja et al. [6], groundwater storage change is associated with vegetation cover, because groundwater supports the growth of natural vegetation and biota. The vegetation water content (VWC) has also been widely used in agriculture and forestry studies. Chen et al. [7] reported that the change of vegetation and VWC have a positive linear relationship according to moderate-resolution imaging spectroradiometer (MODIS) data and ground-measured VWC data. The increase in vegetation will increase soil infiltration, and the high infiltration rate of water into the soil will improve soil water availability [8]. Zheng [9] found that the change of vegetation can also affect soil erosion, using field survey data from Loess Plateau in 2005. Puigdefábregas [10] found that vegetation affects the sediment fluxes in dry lands using field observations and numerical simulation experiments.
The amount of plants on the Earth has increased over the last two decades, and China and India lead the increase in greening on land according to data from the National Aeronautics and Space Administration (NASA) Earth Observatory. This was first reported by Chen et al. [11] based on MODIS measurements. This greening could be expressed as the increase of the leaf area index (LAI). We found that all these mass changes mentioned above related to vegetation changes can be detected by NASA's Gravity Recovery and Climate Experiment (GRACE) gravity satellite.
The GRACE satellite mission has been accurately mapping the variations in Earth's gravity field since its launch in March 2002. The time-varying gravity field solved by GRACE can be used for many applications, for example, in increasing our understanding of climate change and related natural phenomena [12]. The gravity variations studied by GRACE include the mass changes of terrestrial water reserves at large scales (global, regional, and watershed scales) as follows [13,14]: Co-seismic gravitation field variation of large subduction earthquakes [15,16]; global high mountain glacier and polar ice sheet mass change [17][18][19]; global mean sea level budget and sea level change caused by deep ocean warming [20,21]; and sediment accumulation rate change in coastal areas [22]. However, GRACE has not been used to study vegetation mass change.
In the present study, we found that the greening pattern in Southern China is consistent with the mass trend pattern detected by GRACE in Southern China, and we confirmed the greening by GRACE. We analyzed the contribution of each factor and the intra-annual change relationship between vegetation and mass in Southern China using multi-source data.

Leaf Area Index
In this study, we used the MCD15A2H version 6 /level 4 product [23] and MOD17A3 V055 [24] provided by the NASA Land Processes Distributed Active Center Distribution Server hosted at the U.S. Geological Survey Earth Resources Observation and Science Center (https://e4ftl01.cr.usgs.gov/MOTA/). The MCD15A2H product includes the LAI band from the combination of Terra MODIS and Aqua MODIS with a spatiotemporal resolution of 500 m and 8 d. The MODIS LAI products use a sinusoidal projection to grid the data and divide the world into areas measuring approximately 10 • by 10 • . We used the MODIS Reprojection Tool to mosaic the grid areas together and resampled to 0.005 • and converted to geographical projections with the WGS84 datum; we then converted them to GeoTIFF files.

GRACE Data
The GRACE satellite affords researchers an important geodetic method by which to detect the gravity changes on the Earth. It has provided monthly global gravity solutions for more than 10 years since it launched in 2002. We used the monthly Release 06 solutions (spherical harmonic coefficients expanded up to degree 60) provided by the Center for Space Research of the University of Texas at Austin (http://icgem.gfz-potsdam.de/series/01_GRACE_monthly/CSR%20Release%2006). The DDK4 filter was used to process the GRACE coefficients [25]. We added back the degree-1 coefficients according to GRACE Technical Note 13, which were computed based on the data of Sun et al. [26] and replaced the C20s coefficients calculated by the satellite laser-ranging solutions [27]. The glacial isostatic adjustment (GIA) affects all regions of the world, so it should be considered. The GIA effect is corrected by a three-dimensional Geruo13 model [28].

Precipitation Data
Global precipitation data were obtained from the Global Precipitation Climatology Project (GPCP). GPCP is part of the global energy and water cycle exchange activity operating under the auspices of the world climate research program [29]. The global precipitation product was analyzed by integrating the various satellite datasets of ocean and land and a gauge analysis over land [30]. We used the monthly gridded version 2.3 data from 2003 to 2017, which have 2.5 • by 2.5 • spatial resolutions on global ocean and land (http://gpcp.umd.edu/).

Fitting Method
We used the least-squares method to fit the time series to obtain the trend. The semiannual and annual cycles are also considered in the following equation.
where a is constant; b is the trend; c and d are annual and semiannual amplitudes; x is the time; T 1 and T 2 are the semiannual and annual cycles; ϕ 1 and ϕ 2 are the initial phases, respectively; and ε is the residual error.

Mass Change Trend Related to Greening Detected by GRACE
We calculated the LAI trend for China from January 2003 to August 2016; this is shown in Figure 1a. It can be seen that most parts of Eastern China are greening. According to the research of Chen et al. [11], the greening of Southern China is predominantly a result of tree planting programs. The increase of vegetation will bring many benefits to those areas, such as increased biomass accumulation, soil water, ground water, and VWC, and a decrease of soil erosion and sediment. All of these variations will cause a mass change in those areas, which could be detected by GRACE. The trend of mass change from January 2003 to August 2016 in the study area was calculated using GRACE data, as shown in Figure 1b. There are obvious increasing mass signals in Southern China. The statistical significance analysis of the trends within our study region were analyzed and are shown in the Section 2 of the supporting information files ( Figure S1). There are also obvious mass decreasing signals caused by excessive depletion of groundwater in north China [13], as shown in Figure 1b. These large decreasing signals affect our research of vegetation variation. Therefore, we selected the area within the blue box (105.5 • E-117 • E, 24 • N-31 • N, with approximately 9.86 × 10 5 km 2 ) in Figure 1 as our study region, which is not affected by the large decreasing signal.

Effects of the Three Gorges Reservoir and Lakes on the Mass Signals
The mass signals detected by GRACE comprise the overall effects of all the mass changes in this region, and it also includes signals of water level changes of Three Gorges Reservoir and lakes; these sources of noise should be removed using other data.
Wang et al. [31] reported that the water storage changes in the Three Gorges Reservoir of China could be detected by GRACE. Therefore, the reservoir effect should be taken into account. In our study, we used the in situ data to remove the contribution of water storage changes in the Three Gorges Reservoir from the original GRACE data. Solving this problem requires the knowledge of the water level change and its pattern of change. As shown in Figure 2b, the water level can be obtained from the China Three Gorges Corporation (https://www.ctg.com.cn/), and we converted those data to a monthly mean water level. The reservoir spans approximately 600 km along the Yangtze River, from the dam to Chongqing City (Figure 2a), and the maximum water level can reach 175 m. Hence, it is difficult to determine the shape and depth of the reservoir. To obtain the water-surface changing pattern with increasing water levels, we used Shuttle Radar Topography Mission (STRM) 90 m digital elevation data to determine the boundary of the water surface; these data are provided by the Consultative Group for International Agricultural Research Consortium for Spatial Information (CGIAR-CSI, http://srtm.csi.cgiar.org/). According to the monthly mean water level and digital elevation data of the reservoir, we determined the extent of the water surface every month. Thus, the monthly mean water level of every grid point (90 by 90 m) within our study period was calculated. We then calculated the trend of water level change at every grid point within the reservoir from January 2003 to August 2016, and expanded them to the spherical harmonic coefficients up to degree 60. We then, smoothed them with the DDK4 filter, which is consistent with the GRACE process, as shown in Figure 2c. Therefore, the signals caused by the Three Gorges Reservoir were removed from the GRACE solutions.

Effects of the Three Gorges Reservoir and Lakes on the Mass Signals
The mass signals detected by GRACE comprise the overall effects of all the mass changes in this region, and it also includes signals of water level changes of Three Gorges Reservoir and lakes; these sources of noise should be removed using other data.
Wang et al. [31] reported that the water storage changes in the Three Gorges Reservoir of China could be detected by GRACE. Therefore, the reservoir effect should be taken into account. In our study, we used the in situ data to remove the contribution of water storage changes in the Three Gorges Reservoir from the original GRACE data. Solving this problem requires the knowledge of the water level change and its pattern of change. As shown in Figure 2b, the water level can be obtained from the China Three Gorges Corporation (https://www.ctg.com.cn/), and we converted those data to a monthly mean water level. The reservoir spans approximately 600 km along the Yangtze River, from the dam to Chongqing City (Figure 2a), and the maximum water level can reach 175 m. Hence, it is difficult to determine the shape and depth of the reservoir. To obtain the water-surface changing pattern with increasing water levels, we used Shuttle Radar Topography Mission (STRM) 90 m digital elevation data to determine the boundary of the water surface; these data are provided by the Consultative Group for International Agricultural Research Consortium for Spatial Information (CGIAR-CSI, http://srtm.csi.cgiar.org/). According to the monthly mean water level and digital elevation data of the reservoir, we determined the extent of the water surface every month. Thus, the monthly mean water level of every grid point (90 by 90 m) within our study period was calculated. We then calculated the trend of water level change at every grid point within the reservoir from January 2003 to August 2016, and expanded them to the spherical harmonic coefficients up to degree 60. We then, smoothed them with the DDK4 filter, which is consistent with the GRACE process, as shown in Figure 2c. Therefore, the signals caused by the Three Gorges Reservoir were removed from the GRACE solutions. According to the study of Yi et al. [32], the water levels of Dongting Lake and Poyang Lake (see Figure S2) increased in our study period, but the data collection period was from 2003 to 2011, which was shorter than our study period. They found that the levels of Dongting Lake and Poyang Lake are increasing at a rate of 0.02 ± 0.07 and 0.02 ± 0.08 m/yr, respectively, and we assumed that the increased According to the study of Yi et al. [32], the water levels of Dongting Lake and Poyang Lake (see Figure S2) increased in our study period, but the data collection period was from 2003 to 2011, which was shorter than our study period. They found that the levels of Dongting Lake and Poyang Lake are increasing at a rate of 0.02 ± 0.07 and 0.02 ± 0.08 m/yr, respectively, and we assumed that the increased rates of both lakes were 0.02 ± 0.08 m/yr from January 2003 to August 2016; this signal was also removed from the GRACE solutions. When this signal was expanded to the spherical harmonic coefficients up to degree 60 and smoothed with DDK4 filter, the contribution of the two lakes is pretty small. As the influence of these lakes is much smaller than that of Three Gorges Reservoir, we do not show it here. The mean mass change rate of Three Gorges Reservoir and the two lakes is 1.2 ± 0.2 mm/yr in our study region. The mass change pattern shown in Figure 2d with the effects of the reservoir and lakes removed.

Contribution of Each Factor
After removing the noise signals, the mass change in our study region was related to the vegetation variations. The mean mass change rate in our study region was 6.7 ± 0.8 mm/yr from January 2003 to August 2016, as detected by GRACE, which includes the changes of vegetation alone, soil water, and groundwater. The uncertainty of GRACE includes the measurement noise and trend fitting error within a 95% confidence interval. The measurement error was calculated using the method proposed by Wahr et al. [33], and the error of region mean mass trend was obtained by propagating the stokes coefficient errors and trends fitting error.
The annual mass change of carbon in plants was represented by net primary productivity (NPP, the annual accumulation of carbon) which is related to the LAI [34]. Schiesinger and Bernhardt [35] reported that the carbon content of plants is typically approximately 45%-50% (we use the mean value 47.5% ± 2.5%) of the plant tissue (dry). However, the plants contain a lot of water. According to the research of Kramer [36], the water content of herbaceous plants is approximately 80%-90% of the fresh weight, while it is over 50% for woody plants. In our study region, most plants are trees according to the supplementary material of Chen et al. [11], as shown in Figure S3. Thus, the water content of vegetation should be greater than 50%. We used the MOD17A3 data to calculate the mean annual NPP during 2003-2014. As the interannual variation of NPP is very small within our study period, we only show the pattern of NPP in 2014 ( Figure 3). The unit of NPP is kg/m 2 , so we can convert it to equivalent water height using the water density of 1000 kg/m 3 (NPP/1000, unit: m). Then, we expanded it to the spherical harmonic coefficients up to degree 60 and smoothed with DDK4 filter. The mean NPP in our study region is approximately 0.6 mm/yr. Therefore, the total contribution of vegetation ( NPP carbon content × water content = 0.6 (47.5±2.5%) × (≥50%) ) will be larger than 2.5 ± 0.2 mm/yr. Remote Sens. 2020, 12, x FOR PEER REVIEW 6 of 11 reported that the carbon content of plants is typically approximately 45%-50% (we use the mean value 47.5% ± 2.5%) of the plant tissue (dry). However, the plants contain a lot of water. According to the research of Kramer [36], the water content of herbaceous plants is approximately 80%-90% of the fresh weight, while it is over 50% for woody plants. In our study region, most plants are trees according to the supplementary material of Chen et al. [11], as shown in Figure S3. Thus, the water content of vegetation should be greater than 50%. We used the MOD17A3 data to calculate the mean annual NPP during 2003-2014. As the interannual variation of NPP is very small within our study period, we only show the pattern of NPP in 2014 ( Figure 3). The unit of NPP is kg/m 2 , so we can convert it to equivalent water height using the water density of 1000 kg/m 3 (NPP/1000, unit: m). Then, we expanded it to the spherical harmonic coefficients up to degree 60 and smoothed with DDK4 filter. The mean NPP in our study region is approximately 0.6 mm/yr. Therefore, the total contribution of vegetation ( = . ( . . %) ( %) ) will be larger than 2.5 ± 0.2 mm/yr. The increase of vegetation will improve the soil water storage capacity. We used the mean value of four versions of the Global Land Data Assimilation System (GLDAS) hydrological model from NASA and the Climate Prediction Center (CPC) model from National Oceanic and Atmospheric Administration (NOAA) to study the trend of soil moisture. The method of calculation is described The increase of vegetation will improve the soil water storage capacity. We used the mean value of four versions of the Global Land Data Assimilation System (GLDAS) hydrological model from NASA and the Climate Prediction Center (CPC) model from National Oceanic and Atmospheric Administration (NOAA) to study the trend of soil moisture. The method of calculation is described in the supporting information, and the region mean time series of the five models is shown in Figure S4. The regional mean rate is 2.3 ± 1.0 mm/yr.
We also studied groundwater change using the groundwater stress indicators dataset [37] provided by Goethe University, Frankfurt (Germany); the data introduction and data processing methods are described in the supporting information file. The data period was 2001-2010, but we can see that the linear change of groundwater is stable (Figure 4 and Figure S5), so we assumed that the region's mean groundwater increase rates remained the same until August 2016. The groundwater change rate in our region is approximately 0.6 ± 0.1 mm/yr.

Challenge of Calculating Mass Change of All Greening Regions
As mentioned above, large mass redistribution on the Earth can be detected by GRACE, since the mass signal detected by GRACE is the combination of all the mass changes, and the signals from adjacent regions affect each other. The gravity signals are complicated in China, as shown in Figure  1b, and the large decrease of groundwater in northern China affects a large part of the region [38]. For instance, the same problem has been encountered in previous studies of Asia high mountain glacial changes, for example, Matsuo and Heki [39] obtained an ice melting rate of −47 Gt/yr, which is overestimated due to the underestimation of Indian groundwater. Jacob et al. [40] obtained the rate −4 ± 20 Gt/yr, which is underestimated due to the neglect of the positive signal from the interior of the Tibetan Plateau. Yi and Sun [41] reported that the glacier change of high mountain Asia is −35.0 ± 5.8 Gt/yr, calculated using a space domain inverse method. They considered the signals of northern India and the interior of the Tibetan Plateau around the glacier. Therefore, the signals around the The time series and corresponding linear change rate of each factor are shown in Figure 4. We can see that the vegetation mass change is the main factor in our study region (black line). We summarized the change rates of all factors in Table 1. The total change rate of vegetation alone, soil water, and groundwater is approximately ≥ 5.4 ± 1.0 mm/yr, which is consistent with that of GRACE. After removing the contribution of soil water and groundwater, we found that the contribution of vegetation is approximately 3.8 ± 1.3 mm/yr. The change of vegetation is a relatively large quantity that should not be neglected in Southern China.

Challenge of Calculating Mass Change of All Greening Regions
As mentioned above, large mass redistribution on the Earth can be detected by GRACE, since the mass signal detected by GRACE is the combination of all the mass changes, and the signals from adjacent regions affect each other. The gravity signals are complicated in China, as shown in Figure 1b, and the large decrease of groundwater in northern China affects a large part of the region [38]. For instance, the same problem has been encountered in previous studies of Asia high mountain glacial changes, for example, Matsuo and Heki [39] obtained an ice melting rate of −47 Gt/yr, which is overestimated due to the underestimation of Indian groundwater. Jacob et al. [40] obtained the rate −4 ± 20 Gt/yr, which is underestimated due to the neglect of the positive signal from the interior of the Tibetan Plateau. Yi and Sun [41] reported that the glacier change of high mountain Asia is −35.0 ± 5.8 Gt/yr, calculated using a space domain inverse method. They considered the signals of northern India and the interior of the Tibetan Plateau around the glacier. Therefore, the signals around the target are very important to the research results. The exact geographic position and mass change of groundwater in northern China are not available, so it is difficult to remove them from the GRACE data. It then becomes difficult to calculate the mass change of every greening area.

Relation between Intra-Annual Change of Mass and Vegetation
The LAI data were transferred to the spherical harmonic coefficients up to degree 60, and then converted to grid data that are the same as those of GRACE. The mean value of the same month was calculated as shown in Figure 5, in which the blue line is the mass change and the red line the LAI change. According to Figure 5, it can be seen that the LAI and mass increased from January to July and decreased from August to December. However, the LAI is almost the same in July and August and decreased from September to December. Blanken et al. [42] studied the seasonal water exchange above and within a boreal aspen forest using in situ measurement data. They also obtained the seasonal change of water and LAI, and their results are very similar to those of our study. They used the water balance equation (P = E + D + ∆S, where P is precipitation, E evapotranspiration, D vertical drainage, and ∆S is soil water storage variation) to study water exchange. They found that the value of P-E and soil water storage increase before the LAI reaches its maximum value, and those values start to decrease when the LAI reaches the plateau of maximum value, because evapotranspiration becomes the dominant factor. This is consistent with the mass and LAI change in our study region.

Relation between Intra-Annual Change of Mass and Vegetation
The LAI data were transferred to the spherical harmonic coefficients up to degree 60, and then converted to grid data that are the same as those of GRACE. The mean value of the same month was calculated as shown in Figure 5, in which the blue line is the mass change and the red line the LAI change. According to Figure 5, it can be seen that the LAI and mass increased from January to July and decreased from August to December. However, the LAI is almost the same in July and August and decreased from September to December. Blanken et al. [42] studied the seasonal water exchange above and within a boreal aspen forest using in situ measurement data. They also obtained the seasonal change of water and LAI, and their results are very similar to those of our study. They used the water balance equation (P = E + D + ∆S, where P is precipitation, E evapotranspiration, D vertical drainage, and ∆S is soil water storage variation) to study water exchange. They found that the value of P-E and soil water storage increase before the LAI reaches its maximum value, and those values start to decrease when the LAI reaches the plateau of maximum value, because evapotranspiration becomes the dominant factor. This is consistent with the mass and LAI change in our study region. The large seasonal mass change must be related to the rainfall and the seasonal fluctuation of vegetation water content. The mean rainfall of every month in our study region is shown in Figure 5. From January to April, the region rainfall increased much more than the region's mean mass change. This means that most part of the rainwater ran away due to the low soil water storage capacity, which may be related to the low LAI. As the LAI increases, the water storage capacity increases, and the mass increases. When the rainfall starts to decrease, the water storage capacity still increases with LAI increase in June. Thus, we speculate that vegetation has an important effect on the region's water The large seasonal mass change must be related to the rainfall and the seasonal fluctuation of vegetation water content. The mean rainfall of every month in our study region is shown in Figure 5. From January to April, the region rainfall increased much more than the region's mean mass change.
This means that most part of the rainwater ran away due to the low soil water storage capacity, which may be related to the low LAI. As the LAI increases, the water storage capacity increases, and the mass increases. When the rainfall starts to decrease, the water storage capacity still increases with LAI increase in June. Thus, we speculate that vegetation has an important effect on the region's water storage capacity, and, moreover, that it can be detected by GRACE. According to the research of Nemani and Running [43], the LAI and soil available water capacity are positively correlated, and an LAI change of 1 m 2 /m 2 corresponds to a change of approximately 100 mm of soil available water capacity, which is consistent with our study.

Conclusions
In this study, we investigated the detectability of mass change related to variation of vegetation using GRACE data. We used the MODIS LAI dataset to represent vegetation change. The change trends of LAI from January 2003 to August 2016 were large in Southern China. The greening trends in Southern China could be detected by GRACE after removing the effects of the Three Gorges Reservoir and nearby lakes. The mass change in our study region was 6.7 ± 0.8 mm/yr from January 2003 to August 2016 detected by GRACE, and the sum of each factor calculated by multi-source data is consistent with GRACE. Therefore, the greening of Southern China was confirmed by GRACE. The mass change of vegetation is approximately 3.8 ± 1.3 mm/yr calculated by GRACE, and thus, non-negligible in the study of other scientific problems. We also analyzed the intra-annual change of vegetation and mass, and we determined that there is a strong positive correlation between the intra-annual change of mass and vegetation. Furthermore, the change of vegetation could regulate the water storage capacity of land, a feature was also detected by GRACE.
Supplementary Materials: The following are available online at http://www.mdpi.com/2072-4292/12/2/328/s1, Figure S1: Distribution of p-values across our study region, Figure S2: Positions and water levels of two lakes in southeastern China, Figure S3: Map of vegetation classes aggregated from 2007 MCD12C1 product, Figure S4: Time series of regional mean soil water change, Figure S5: Time series of regional mean groundwater change from 2003 to 2010.