Quantitative analysis of the impacts of terrestrial environmental factors on precipitation variation over the Beibu Gulf Economic Zone in Coastal Southwest China

Taking the Guangxi Beibu Gulf Economic Zone as the study area, this paper utilizes the geographical detector model to quantify the feedback effects from the terrestrial environment on precipitation variation from 1985 to 2010 with a comprehensive consideration of natural factors (forest coverage rate, vegetation type, terrain, terrestrial ecosystem types, land use and land cover change) and social factors (population density, farmland rate, GDP and urbanization rate). First, we found that the precipitation trend rate in the Beibu Gulf Economic Zone is between −47 and 96 mm/10a. Second, forest coverage rate change (FCRC), urbanization rate change (URC), GDP change (GDPC) and population density change (PDC) have a larger contribution to precipitation change through land-surface feedback, which makes them the leading factors. Third, the human element is found to primarily account for the precipitation changes in this region, as humans are the active media linking and enhancing these impact factors. Finally, it can be concluded that the interaction of impact factor pairs has a significant effect compared to the corresponding single factor on precipitation changes. The geographical detector model offers an analytical framework to reveal the terrestrial factors affecting the precipitation change, which gives direction for future work on regional climate modeling and analyses.


Results and Interpretation
Precipitation change in the Beibu Gulf Economic Zone. The precipitation trend of each meteorological observation (88 total meteorological observatories in Guangxi Province) was calculated from annual precipitation for the period of 1985~2010 using Equation (1). Linear trends indicate that 76% of the meteorological stations show a positive trend in annual precipitation during 1985-2010, and notably, five of them are statistically significant at the 90% confidence level (see Supplementary Table S1, Fig. S1). The other 24% of the meteorological stations show a negative trend (see Supplementary Table S1). A station-by-station analysis was performed and mapped using ArcGIS 10.1 with the Empirical Bayesian Kriging interpolation method in order to explore spatial patterns of precipitation changes in Guangxi Province. Then, the precipitation changes in the Guangxi Beibu Gulf Economic Zone were clipped and are shown in Fig. 2. This might be better than using direct interpolation of fewer meteorological observatories from the Guangxi Beibu Gulf Economic Zone, especially in the border area, because of the regional characteristics of precipitation.  Figure 2 shows that the precipitation trend rate of the Beibu Gulf Economic Zone is between − 47 and 96 mm/10a. The character of the spatial precipitation change is primarily in the northwest-southeast direction. The low value zones are located in the southwest Beibu Gulf Economic Zone, while the high values are located in the northeast. The figure also shows that the precipitation changes in the middle zone are relatively lower than those in the neighboring south and north areas. A relative increasing precipitation trend from south to north is observed as a whole. We also found two zero lines of precipitation change in the middle and southwest zone of the Beibu Gulf Economic Zone.
According to Li & Su's research, the Mann-Kendall mutability test found that precipitation in Guangxi Province had sudden changes in the years of 1984 and 1994. Specifically, from 1984 to 1994, Guangxi Province had less rain, while beginning in 1994, Guangxi Province entered into a relatively pluvial period 44 . Existing research also shows that there was a positive trend center of precipitation in northwest Guangxi. Therefore, in general, the linear trend is increasing.
The feedback of terrestrial environmental factors to precipitation change. The leading factors of precipitation change. The factor detector ranked the terrestrial environment layers by their influences (P D,H values) on precipitation change in the following order for the study area (Table 1): Among the terrestrial environmental factors, the P D,H value of FCRC is the maximum. Obviously, there is a large break in sorted P D,H values between PDC and VT. The P D,H values of FCRC, URC, GDPC, FRC and PDC are in a group with high values and small differences, while the rest of the factors belong to another group with lower values. In a general sense, if the P D,H value of a factor is larger than 0.2(20%), then the factor can be regarded as a leading factor 39 that strongly explains the spatial pattern. Therefore, FCRC, URC, GDPC, FRC and PDC are potential leading factors that may exert the largest impact on precipitation change (spatial pattern) in this study area. In contrast, the values of VT, DEM, GT, GRD, ASP, LUCC and TET are comparatively small at less than 0.1(10%), which likely reflects their smaller contributions to the precipitation trend (spatial pattern).
The ecological detector (Table 2) shows the differences of the P D,H values. Among the five potential leading factors (FCRC, URC, GDPC, FRC and PDC), approximately 80% of them (FCRC, URC, GDPC and PDC) are not statistically significant with each other, whereas statistically significant differences between FRC and other potential leading factors (URC and GDPC) were found. That is, URC and GDPC have a larger significant effect on the precipitation change than FRC. With the factor detector and the ecological detector, we concluded that FCRC, URC, GDPC and PDC are leading factors, and FRC was eliminated from the potential leading factors. Therefore, FCRC, URC, GDPC and PDC have the largest contribution to the precipitation change, whereas the remaining factors have a relatively weak influence.  In view of the above considerations, we also found that the power of social factors is much larger than that of natural factors in changing precipitation in the study area. This can be seen from the leading factors (FCRC, URC, GDPC and PDC) and whole sorting of the P D,H values, which means that people are likely to be a very powerful factor in changing the terrestrial environment to influence local precipitation on a regional scale.
On the surface, FCRC is the first leading natural factor. However, the increase of forest coverage rate change in the study area is predominantly due to ecological construction for tourism and sustainable development in recent years. The forest coverage rate change of this study area is higher than the rate of other places in China and occurs in a sustainable growth manner. The Shiwan Dashan National Forest Park and the Daming Mountain National Natural Reserve are located in this research area. Some research has already highlighted that the significant role of precipitation may increase or decrease alongside afforestation 4-6 or deforestation.
The Guangxi Beibu Gulf Economic Zone is the first international regional economic cooperation zone in China. According to the statistics, the population of this region was 15.81 million in 1985 and rose by 43 percent to 22.68 million in 2010. The GDP of this region increased from 8.3 billion Yuan (RMB) to 412.1 trillion Yuan (RMB) between 1985 and 2010, which is a huge growth of 49.37 times the initial GDP. With the sustainable growth of population and the increasing development of the economy, the GDP, the population density, urbanization, construction activities, energy consumption and greenhouse gas emissions have experienced a relatively rapid growth in the Beibu Gulf Economic Zone, and among them, the GDP growth rate is the largest.
The urbanization levels of both Nanning City and Beihai City are over 55%. In contrast, Guangxi has a relatively extensive development model for their economy, with an energy consumption per unit of GDP of 1.036 tons of standard coal per ten thousand Yuan (2010), which is 1.28 times the national average of 0.81 tons of standard coal per ten thousand Yuan. A large amount of energy consumption emits a large amount of greenhouse gas, such as carbon dioxide, which is the main source 45 of carbon emissions. This greatly influences the climate of this region and possibly the climate on a larger regional scale. Urbanization is a comprehensive process, which will influence a city's precipitation, temperature, humidity, visibility and wind, forming a special local meteorological environment and causing material climate changes. The GDP may actually be viewed as a comprehensive result of many human activities. The increase of population density leads to an increase of artificial thermal discharge, directly influencing the change of surficial sensible heat flux, which will influence precipitation significantly 46 . In addition, approximately 45% of China's farmland is irrigated farmland 47 , whereas Guangxi Province has a higher percentage. The P D,H value of FRC on precipitation change is 35.3% (much higher than 20%). This occurs mainly because the heavy irrigation of farmland affects the distribution of surface net radiation between latent heat flux and sensible heat flux change (latent heat flux increases, but sensible heat flux decreases), and farm irrigation has a cooling effect on the earth's surface; at the same time, the increase of soil humidity enhances transpiration and further increases the moisture content in the atmosphere and the unstable energy of latent heat, leading to an increase of convective precipitation 48 and producing a marked effect on the region's precipitation. This finding is supported by other cases. Irrigation over the Ogallala Aquifer in the central United States increased dramatically over the 20th century and has enhanced regional precipitation 49 . The precipitation increase in the Texas Panhandle from 1952 to 1980 was obviously due to the increase in the irrigation area 25 . On the other hand, the amount of precipitation in central and southern India decreased due to a lower surface temperature over the irrigated areas of India in July 26 .
The Beibu Gulf Economic Zone is a relatively small area, which is on a small scale compared to the majority of research on precipitation. In the region, the DEM, geomorphic type, slope aspect, gradient, ecosystem and vegetation form are similar or experience less change, so they probably have a weak effect on precipitation change.
The effect of the interaction of terrestrial environmental factors on precipitation change. The interaction detector was used to check whether or not two factors work independently. The joint impacts of two factors measured by the P D,H values are shown in Table 3 and Table S2 and can be compared with their separate impacts.
It must be noted from Table S2 that the P D,H values of 22 interactive pairs are greater than that of the primary leading factor (FCRC). The max P D,H value comes from interaction of FRC with GDPC (FRC ∩ GDPC = 84.4%). Specifically, all the interactive effects between FCRC and the rest of the factors (FCRC ∩ GDPC = 83.1%,  , and all of them are larger than 70%. We thought that FCRC is also a social factor because FCRC is mainly due to human ecological construction. Therefore, these factors are all social factors, and it clearly implies that humans are the most important aspect in changing precipitation (similar to the analysis of leading factors) in this region via economic activities such as urban construction, afforestation, changing and developing hillside fields, irrigation and plantation. Under the high pressure of growing population and development, humans are the best medium compared to other natural factors to change and affect the spatial distribution of other factors according to their purposes, and with the development of science and technology, this situation is amplified. For example, large-scale afforestation in the northern mid-latitudes warms the Northern Hemisphere and alters global circulation patterns to redistribute the anomalous energy absorbed in the northern hemisphere, which results in a precipitation decrease over parts of the Amazon basin and an increase over the Sahel and Sahara regions in Africa 50 .
In addition to the above mentioned, we also noted that interactions between social factors and natural factors have two types: nonlinear enhancement and bienhancement (Table 3). Each type indicates that the factors bienhance or nonlinearly enhance each other. As shown in Table 3, the interactions between social factors and natural factors have predominantly strong, nonlinear synergies. For example, the interactions of PDC and DEM (PD ∩ DEM = 40.3% > 34.7% = PDC (27.4%) + DEM (7.3%)) are larger than the P D,H value sum of PDC and DEM; therefore, the interaction between PDC and DEM has a larger impact on precipitation changes. This is likely due to the city and farmland expansion toward a relatively bad condition of DEM that changes the underlying surface conditions. It also indicates that social factors and natural factors have synergies and can enhance each other's effect on precipitation change.
In conclusion, social factors have a larger impact on the precipitation change compared to natural factors. Partial natural factors have a relatively small impact on precipitation change but show a strong synergy with the interaction of other factors. The feedback of terrestrial environmental factors on precipitation change mainly arises from interactions of impact factors and interactive pairs of impact factors, which have a larger influence on precipitation change than the single factor does through the feedback. Interactions between factors play an important role in the precipitation change in this region.

Regional analysis of the leading impact range (type) of leading factors on precipitation change.
The risk detector shows that the average precipitation change in the different FCRC zones (from I to VI) are − 0.66 mm/10a, − 9.98 mm/10a, 53.18 mm/10a, 6.88 mm/10a, 23.48 mm/10a and 11.66 mm/10a, respectively, and they are significantly different. It also implies that precipitation will increase or decrease with the increase or decrease of forest coverage. However, higher precipitation change is not consistent with a larger FCRC zone, and precipitation change fluctuates with FCRC values. A similar analysis of other terrestrial environmental factors can be conducted using the risk detector. The small and continued growth of annual urbanization rates will lead to a large increase in annual precipitation. The main impact ranges of FCRC, URC, GDPC and PDC tend to be located at the relatively low-middle value zones. We selected the largest types (ranges) of each leading factor as the main impact types (ranges) by sorting the average precipitation change. The main impact types (range) are tabulated in Table 4 and mapped in Fig. 3.
As shown in Fig. 3, the leading impact type or range of each leading factor on precipitation change is predominantly located in the northeast-southeast of the Beibu Gulf Economic Zone. This indicates that the largest precipitation change is in the northeast-southeast of the Beibu Gulf Economic Zone, and the range of the precipitation trend rate is between 39 and 96 mm/ 10a (Fig. 3). Therefore, the main distribution areas of the main impact range (type) of the leading factors on precipitation revealed by the results of the risk detector are consistent with the distribution of the relatively large area of precipitation change trend rate calculated by the linear regression model. This illustrates the flexibility of applying the geographical detector model to obtaining initial detection results of the precipitation change mechanism. Figure 3 shows that the precipitation change for the county of Rongxian is strongly controlled by FCRC, URC and PDC. According to the interaction detector, we also found that the P D,H values of FCRC ∩ URC (74.1%), FCRC ∩ PDC (75.0%), and URC ∩ PDC (71.0%) are very high and enhance each other to increase precipitation change, which emphasizes directions for future work. In conclusion, the largest precipitation change is present in the northeast-southeast region of the Beibu Gulf Economic Zone and is predominantly influenced by the interactions of factors such as FCRC, URC, GDPC and PDC.

Conclusions and Discussion
The causes of precipitation changes are very complicated due to the interaction of the land surface with the atmosphere. In addition, the research resources, such as shared data, are limited in developing countries, creating a high demand for useful detecting and/or analyzing tools. In this study, we used geographical detectors to verify the effects of some of the natural and social factors on precipitation change at a regional scale. We believe that this program is unique because it extracts the interrelationships between precipitation change and terrestrial environmental factors using the correspondence of their spatial distribution and, most importantly, because it is easily implemented. The feedback of terrestrial environment to precipitation changes can be partially explained by forest cover, urbanization, terrain, irrigation and other single factors. Typically, the comprehensive consequences are the result of interactions of multiple factors. In this study, we found the following:  Our research suggests that the geographical detector offers a quantitative and objective analytical framework that could be used to find the essence of many geosciences phenomena. There are still several aspects for future study. First, spatial scale transformation is an important aspect of geographical detectors. Transforming the administrative regions into the same grid cells might be subjective, as the grid size can have different values. We also found that discretization methods to classify continuous variables into several categories might affect the results because these methods do not currently have standardized rules. Second, due to the limitation of range and data accessibility in this study area, quantitative analysis was not conducted overall based on impact factors in this study. Third, the main impact ranges of leading factors (FCRC, URC, GDPC and PDC) fluctuated with precipitation change, and the largest precipitation change is typically only consistent with the smallest range of URC. In the future, threshold values of the main impact ranges can be overcome by collaboration with climate models. This is likely a better way to integrate geographical detectors with traditional meteorology methods to discover the precipitation change mechanism.
Despite some limitations, we still believe that this study will be meaningful. The geographical detectors are statistical and are not a causality tool; however, they can distinguish high potential impact factors and leading factor ranges to emphasize the next step in research. The results from this study can help researchers to understand the spatial pattern of precipitation change with impact factors and provide clues for further studies by integrating traditional observation, simulation, contrast testing, etc.

Materials and Methods
Research methods. Trend rate. Tests for trend detection of the climatic element in a time series can be classified as parametric and non-parametric methods (e.g., the Mann-Kendall test). The linear regression method is a very simple and common parametric method 51 , and the trend rate method generally adopts the unitary linear regression model, that is: where y represents a climatic element or other sequence (e.g., precipitation); x represents a yearly time series (from 1985 to 2010); and b represents a linear trend term, the value of which is a linear trend rate, in mm/10a.
Geographic detector model. Geographical detectors are composed of the factor detector, ecological detector, risk detector and interaction detector 34,43 . Factors significantly affecting precipitation change can be selected as the leading factors through analysis using the factor detector and ecological detector models; the risk detector can further analyze leading impact types or scopes (confidence level of 95%) of impact factors that significantly affect precipitation change; and the interaction detector can analyze the interaction among various factors. The core concept of the factor detector is as follows: there is certain differentiation of the factors affecting the development of geographical phenomenon in space. If a certain factor has a remarkable consistency with the change of that geographical phenomenon in space, then the factor will have a definite determinant power on the occurrence and development of a geographical phenomenon 34 , measured by the size of the power determinant value (P D,H ). Details of the geographical detector can be found in the original paper 34 . Here, in our research context, the calculation model for detecting impact factors of precipitation change in the Beibu Gulf Economic Zone is reviewed as follows: We assume that precipitation change would present a spatial distribution similar to that of an impact factor if the impact factor leads to the change of precipitation (see Supplementary Fig. S2). All impact factors are quantified by these power values as follows: In equation (2), D represents an impact factor layer (e.g., DEM or slope) that must already be categorized (e.g., DEM values can be categorized into eight categories); m is the number of zones (categories) of the factor D (D = {D 1 ,D 2 ,D 3 , … , D m }); H represents the precipitation trend rate; P D,H represents the power of determinant D on H; n and σ H 2 represent the number of total samples and the global variance of H over the entire study area, respectively; n D i , and σ H The ecological detector compares which suspected impact factor (e.g., C factor) determinant is more significant than the other (e.g., D factor) in causing precipitation change in the study area. This is measured using the F-test: ( 1 ) In equation (3), F is the test value of F, n C,p and n D , p denote the number of samples of impact factors C and D in sample unit p, respectively, and σ C,m 2 and σ D,m 2 are dispersion variances of impact factors C and D, respectively. The null hypothesis is H 0 : σ C,m 2 = σ D,m 2 . If H 0 is rejected conditioned on a significant level α (usually 5%), we conclude that the impact factor C is more significant than the impact factor D in affecting precipitation change.
Scientific RepoRts | 7:44412 | DOI: 10.1038/srep44412 (d) Geomorphic type derived from a 1: 1,000,000 geomorphic map at the spatial resolution of 1000 * 1000 m; (e) Data of the terrestrial ecosystem type were derived from spatial distribution data of the Chinese terrestrial ecosystem types at the spatial resolution of 1000 * 1000 m; (f) Forest coverage rate change was calculated through Equation (1), and then, the trend of the forest coverage rate of each county (b value) was mapped and classified into 6 types through NB; (g) Data of the vegetation types were derived from a 1: 1,000,000 vegetation map at the spatial resolution of 1000 * 1000 m. Different types or ranges of an impact factor have different influences on precipitation change. The risk detector compares the differences through the t-test. The computational formula is as follows: In equation (4), t ij is the test value of t; R i and R j are average values of the precipitation tendency rate over property i and property j of the impact factor R; σ i 2 and σ j 2 are the variances of the precipitation tendency rate from property i and property j, respectively; and n i and n j are the sample sizes of the two properties.
The interaction detector shows that when the two different factors of x and y are combined, they either weaken or enhance each another or they are independent in changing precipitation, determined by comparing P D,H (x ∩ y) with the values of P D,H (x) and P D,H (y), where the symbol ' ∩ ' denotes the intersection between the x layer and y layer. If P D,H (x ∩ y) < min (P D,H (x), P D,H (y)), the variables nonlinearly weaken each other; if min (P D,H (x), P D,H (y)) < P D,H (x ∩ y) < max (P D,H (x), P D,H (y)), the variables uniweaken each other; if P D,H (x ∩ y) > max (P D,H (x), P D,H (y)), the variables bienhance each other; and if P D,H (x ∩ y) > P D,H (x) + P D,H (y), the variables nonlinearly enhance each other. If P D,H (x ∩ y) = P D,H (x) + P D,H (y), then the variables are independent of each other.
Based on the precipitation trend rate in the Beibu Gulf Economic Zone from 1985 to 2010, geographical detectors were utilized to explore the impact and indication effect of terrestrial environmental factors on the precipitation change through climate feedbacks.
Technical process. Modeling of the geographical detector mainly involves the following steps: first, determination of the optimal classification method for the factor data; second, determination of the impact of factors on the precipitation change; and third, determination of the leading role of factors in the precipitation change. Regarding the technical process in detail, please see Supplementary Fig. S3.
Data sources and processing. Precipitation data. The selected observation data are the mean annual precipitation of 29 meteorological stations from 1985 to 2010 in the Beibu Gulf Economic Zone (Fig. 1). The above data were derived from the China meteorological data network (http://data.cma.cn).
Potential natural factors. According to the main impact factors of precipitation change discussed in the introduction, almost all of the environmental factors, except for climate type, were considered as main potential natural factors, such as geomorphic type, the types of terrestrial ecosystem, vegetation type, elevation, gradient, aspect and forest coverage rate, to reveal the feedback. Based on the results of China's ecological geographic division, the entire study area belongs to the climate type of the south subtropical-humid region; thus, the climate type factor can be ruled out here.
During the study period , the elevation, gradient and aspect of the Guangxi Beibu Gulf Economic Zone remain relatively stable, so SRTM DEM was used and also to produce the gradient and aspect. Similarly, the changes of geomorphic type, the type of terrestrial ecosystems and vegetation type were relatively small and fragmented, so we selected a middle year (around 2000) of these datasets to represent the entire study period. The datasets above were provided by the Data Center for Resources and Environmental Sciences, Chinese Academy of Sciences (RESDC) (http://www.resdc.cn). Annual forest coverage rates were collected from the Guangxi Forestry Yearbooks  and the Guangxi Statistical Yearbooks. To ensure the continuity of the dataset in time, regression analysis was used to fix missing data.
Data sorting and pretreatment were conducted in ArcGIS10.1. Based on the input requirements of the geographical detector model (http://www.sssampling.org/Excel-GeoDetector/), projection was unified to the projection coordinate system of Krasovsky-1940-Albers, and raster data were reclassified as 6 to 8 grades 27,36,37 and then converted to the vector data type. ArcGIS provided some discrete classification methods, such as the Equal Interval Method (EI), Quantile Value Method (QV), Natural Break Method (NB) and Geometrical Interval Method (GI), to reclassify the raster data. Different classification methods result in different P D,H values for the classified factor. The highest P D,H value result will indicate that this impact factor classification, using the discrete method, can be more representative as the classification of a geographical phenomenon, thus better revealing spatial distribution laws of the geographical phenomenon 38 . Natural factors were processed and classification methods were selected after many experiments (Fig. 4).
Potential social factors. Population density, GDP, farmland rate, urbanization rate and land use were selected as potential social factors that likely caused regional precipitation change because of changes in them, as described in the introduction. The population density, GDP and urbanization rate were derived from the Guangxi Statistical Yearbooks (1986~1991, 1993~1999 and 2001~2010), while the farmland rate was derived from the Guangxi Rural Statistical Yearbooks (1985~2010), and regression analysis methods were used to fill the entire 26-year period . Land use data (1980s and 2010) were collected from the Data Center for Resources and Environmental Sciences, Chinese Academy of Sciences (RESDC) (http://www.resdc.cn). We used the trend rates of these factors, derived from Equation (1), to express change because these factors have changed greatly over the 26-yearperiod. Similarly, the social factors adopted the same processing method as the potential natural factors (Fig. 5).