ORIGINAL RESEARCH published: 23 January 2020 doi: 10.3389/fdata.2019.00051 Spatial Pattern and Environmental Drivers of Acid Phosphatase Activity in Europe Yan Sun 1*, Daniel S. Goll 1,2 , Philippe Ciais 1 , Shushi Peng 3 , Olga Margalef 4,5 , Dolores Asensio 4,5 , Jordi Sardans 4,5 and Josep Peñuelas 4,5 1 Laboratoire des Sciences du Climat et de 1’Environnement, CEA-CNRS-UVSQ, Gif sur Yvette, France, 2 Institute of Geography, University of Augsburg, Augsburg, Germany, 3 Sino-French Institute for Earth System Science, College of Urban and Environmental Sciences, Peking University, Beijing, China, 4 CSIC, Global Ecology Unit, Centre de Recerca Ecològica i Aplicacions Forestals, Consejo Superior de Investigaciones Científicas, UAB, Bellaterra, Spain, 5 Centre de Recerca Ecològica i Aplicacions Forestals, Cerdanyola del Vallès, Spain Acid phosphatase produced by plants and microbes plays a fundamental role in the recycling of soil phosphorus (P). A quantification of the spatial variation in potential acid phosphatase activity (AP) on large spatial scales and its drivers can help to reduce the uncertainty in our understanding of bio-availability of soil P. We applied two machine-learning methods (Random forests and back-propagation artificial networks) to simulate the spatial patterns of AP across Europe by scaling up 126 site observations Edited by: of potential AP activity from field samples measured in the laboratory, using 12 Forrest M. Hoffman, environmental drivers as predictors. The back-propagation artificial network (BPN) Oak Ridge National Laboratory (DOE), method explained 58% of AP variability, more than the regression tree model (49%). United States In addition, BPN was able to identify the gradients in AP along three transects in Europe. Reviewed by: Jitendra Kumar, Partial correlation analysis revealed that soil nutrients (total nitrogen, total P, and labile Oak Ridge National Laboratory (DOE), organic P) and climatic controls (annual precipitation, mean annual temperature, and United States Qingyuan Zhang, temperature amplitude) were the dominant factors influencing AP variations in space. Universities Space Research Higher AP occurred in regions with higher mean annual temperature, precipitation and Association (USRA), United States higher soil total nitrogen. Soil TP and Po were non-monotonically correlated with modeled *Correspondence: AP for Europe, indicating diffident strategies of P utilization by biomes in arid and humid Yan Sun ysun@lsce.ipsl.fr area. This study helps to separate the influences of each factor on AP production and to reduce the uncertainty in estimating soil P availability. The BPN model trained with Specialty section: European data, however, could not produce a robust global map of AP due to the lack This article was submitted to Data-driven Climate Sciences, of representative measurements of AP for tropical regions. Filling this data gap will help a section of the journal us to understand the physiological basis of P-use strategies in natural soils. Frontiers in Big Data Keywords: back-propagation artificial network, Europe, phosphorus cycling, partial correlation analysis, Received: 31 May 2019 regression tree, soil acid phosphatase Accepted: 26 December 2019 Published: 23 January 2020 Citation: INTRODUCTION Sun Y, Goll DS, Ciais P, Peng S, Margalef O, Asensio D, Sardans J and Phosphorus (P) is an essential nutrient for all living organisms (White and Hammond, 2008). As Peñuelas J (2020) Spatial Pattern and Environmental Drivers of Acid weathering of minerals and the deposition of atmospheric dust are minor sources of P (Walker Phosphatase Activity in Europe. and Syers, 1976; Vitousek et al., 2010; Wang et al., 2015), the recycling of organic P from litter Front. Big Data 2:51. and soil organic matter is of utter importance for plant growth and microbial activity in terrestrial doi: 10.3389/fdata.2019.00051 ecosystems. In P-poor ecosystems, limited P recycling may dampen the response of plant growth to Frontiers in Big Data | www.frontiersin.org 1 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity elevated CO2 concentration (Ellsworth et al., 2017). Yang et al. Keenan et al., 2012; Were et al., 2015), but have not been used for (2014) reported in a modeling study that the effect of elevated upscaling the spatial patterns of AP. CO2 on plant productivity in the Amazon Basin critically We used two ML methods in combination with gridded dependent on assumptions regarding the P-recycling efficiency fields of environmental factors to upscale site data of potential within soils, which was strongly related to the parameterization AP (Margalef et al., 2017) to gridded AP fields for continental of phosphatase production in their model. Europe. Then we identify the main drivers behind the spatial The rate at which ecosystems can recycle P from litter and variation of upscaled AP in Europe. Finally, we used the best ML soil organic matter is poorly quantified by observation (Gill model trained by European data in a first attempt to produce and Finzi, 2016). Soil phosphatases secreted by fungi, bacteria, a global map of AP, and this map is cross-validated using non- and plant roots play an important but poorly quantified role in European data. transforming complex and unavailable forms of organic P into assimilable phosphate (Caldwell, 2005). Potential phosphatase METHODS activity in soils, which can be measured in the lab from soil samples, is an indicator of the capacity of enzyme communities AP measurements compiled by Margalef et al. (2017) are to cleave organic molecules containing P (Krämer and Green, combined with gridded fields of 12 environmental factors to 2000), and serves as a surrogate for the lacking measurements of model spatial pattern of annual AP across Europe and the P mineralization in the soil. globe on a spatial resolution of 10 km. This was done by using Potential phosphatase activity measured in the laboratory back-propagation artificial network (BPNs) and Regression Trees under optimal conditions (optimal temperature, well-mixed soil, (RT). Figure 1 shows the main procedure of model calibration no water limitation; Eivazi and Tabatabai, 1977) provides an and evaluation. upper limit of their actual activity in a soil (Nannipieri et al., 2011; Margalef et al., 2017), which cannot be directly measured. Datasets Acid phosphatases (AP) are more widespread than alkaline We used 296 measurements of AP from 139 published studies phosphatases at soil pH values representative of most natural (Margalef et al., 2017). This dataset contained 54 tropical sites and soils (Margalef et al., 2017), which justifies our focus on AP in 242 ex-tropical sites (155 in Europe). Location, soil pH and soil this study. nutrient contents (total soil C, soil organic C, total nitrogen, and Many experiments have investigated the responses of potential total P) for each site were obtained from the original publications. AP activity to fertilization (Marklein and Houlton, 2012; Maistry We separated the data set into four subsets as shown in et al., 2015), temperature changes (Sardans et al., 2006), and Figure 1: (1) data used to train the ML models (Data set A), (2) water availability (Sardans and Peñuelas, 2005; Sardans et al., data used for cross-validating the spatial gradients of AP modeled 2006; Zhou et al., 2013) or to other disturbances under controlled for Europe (Data set B), (3) data used to evaluate AP patterns in conditions (Sinsabaugh et al., 2008). Gradients of AP along temperate regions outside Europe (Data set C), and (4) data used transects have been measured in few regions (Brockett et al., to evaluate AP patterns in tropical regions (Data set D). Dataset A 2012; Huang et al., 2013; Kitayama, 2013). They found that contains 126 European sites with complete information for the all warming, increasing soil water and nitrogen availability can the predictors listed below. Data set B contains all 155 European enhance AP activity. However, these studies were limited to site sites. Data set C and Data set D contained 87 temperate sites or small region scale and only considering a subset of potential outside Europe and only 54 tropical sites, respectively. environmental factors. Recently, Margalef et al. (2017) addressed We selected 12 variables as predictors for the upscaling this gap, by compiling a global data set of phosphatase activity based on findings from manipulation experiments (Hogan and using correlation analysis, regression analysis and structural et al., 2010; Huang et al., 2011; Marklein and Houlton, 2012), equation models (SEMs) to provide insights of the drivers of observation time series (Sardans et al., 2006) and previous phosphatase activity distribution on global scale. However, the regression analyses (Margalef et al., 2017): soil pH (SoilpH), approach used by Margalef et al. (2017) cannot account for soil clay content (Clay), soil organic carbon (OC), soil total non-linear responses to different variables (Ma et al., 2017) and nitrogen (TN), soil total phosphorus (TP) and soil labile organic omitted important variable like soil labile P. phosphorus (Po), net primary productivity (NPP), mean annual The benchmarking of the growing number of land models temperature (MAT), amplitude of yearly temperature (AMP), which include phosphorus cycling is currently hampered by the mean annual precipitation (MAP), soil type (SoilType), and lack of spatial explicit information on AP on regional to global vegetation type (VegType). scale. The limited understanding of the drivers responsible for We extracted predictors from each original publication differences in AP across different ecosystems and climatic and when available. In case information on predictors were not soils conditions further hampers the global efforts of including reported, we extracted the missing data from gridded datasets P cycles in the land surface models (Reed et al., 2015). (Table 1) based on the geographical coordinates of each Machine learning (ML) is a family of approaches which measurement site (Shangguan et al., 2014; Ballabio et al., 2016, has been increasingly used to identify patterns in complex 2019; Fick and Hijmans, 2017; Hengl et al., 2017; MOD12Q2 ecological data sets and scale up site measurements (e.g., Papale and MOD17). The databases for each predictor variables are and Valentini, 2003; Golubiewski, 2006; Wiesmeier et al., 2011; listed in Table 1. Detailed information about the databases Frontiers in Big Data | www.frontiersin.org 2 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 1 | Flow chart for the calibration and evaluation of the prediction models. TABLE 1 | Predictors used in prediction models and its source. Predictor name Abbreviation Source Spatial References resolution Soil organic carbon OC Soilgrids 250 m Hengl et al., 2017 Soil pH SoilpH LUCAS 500 m Ballabio et al., 2019 Soil total nitrogen TN Soil total phosphorus TP Soil clay content Clay Ballabio et al., 2016 Soil labile organic phosphorus Po Global maps of the soil P contents for different 1 km Yang et al., 2013; Hengl et al., 2017; Sun P forms et al., 2017; Ballabio et al., 2019 TP from LUCAS USDA soil types and ratio of labile inorganic P and labile organic P Net primary productivity NPP MODIS-NPP (MOD17A3; mean value during 1 km Running et al., 2004; Zhao et al., 2005; 2000–2014) Turner et al., 2006 Mean annual temperature MAT WorldClim 1 km Fick and Hijmans, 2017 Yearly temperature amplitude AMP Mean annual precipitation MAP Soil type (categorical variable) SoilType Soilgrids and USDA class 250 m Hengl et al., 2017 Vegetation type (categorical variable) VegType MODIS 1 km Friedl et al., 2010 are provided in Supporting Information Appendix S1. error is achieved between the predicted and observed outputs Details on the gap filling of predictors can be found in (Kelley, 1960). We applied a BPN constituting a four-layer neural Supporting Information Appendix S2. AP, NPP, AMP, network: one input layer, two hidden layers and one output layer. and MAP, which followed a lognormal distribution, Two categorical variables (vegetation type and soil type) were were log-transformed before build the prediction model converted into vector form as BPN are not designed to handle (Supporting Information Appendix S2, Figure S2). categorical values following Mason et al. (2017). For example, the soil type has eight possible categories: (1) Entisol, (2) Inceptisol, Back-Propagation Artificial Networks (3) Aridsol, (4) Mollisol, (5) Alfisol, (6) Spodosol, (7) Ultisol, (8) Back-propagation training algorithm (Kelley, 1960) is the most Oxisol. An Entisol would have a vector value of [1 0 0 0 0 0 0 frequently used neural-network method (Were et al., 2015). 0], whereas a Alfisol would have a vector value of [0 0 0 0 1 0 These algorithms train networks until some targeted minimal 0 0]. The input layer contained 34 nodes (10 numeral predictors Frontiers in Big Data | www.frontiersin.org 3 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity and two categorical predictors; Table 1). The numbers of neurons Extrapolation of AP for Europe for the two hidden layers were 10 and 5 for building the BPN We resampled all of those gridded datasets of predictors (Table 1, model (Supporting Information Appendix S3, Figure S3). Each Table S1) on a 10 km resolution using area-weighted mean layer of the BPN was linked to prior and forward layers by methods for numerical predictors, and by using mode resampling weights that were determined using a gradient-descent learning for categorical predictors (i.e., vegetation types and soil types). algorithm, such as Widrow-Hoff learning rule (Kelley, 1960). We Then we used the 25 selected model ensembles (25 best models × used the resilient back-propagation algorithm to update weights 126 LOO rounds) that produced AP with the highest correlation and biases along the negative of the gradient of the performance with observed AP in combination with gridded maps of the function (Riedmiller and Braun, 1993). predictors at a resolution of 10 km (Figure 2) to scale up point Each BPN runs multiple epochs until either of three observations to continental Europe (36–75◦ N and 20◦ W to convergence criteria is satisfied: (1) the number of training 32◦ E). We first used the median of the 25 estimates for each LOO epochs reach up to 500, (2) the predicting bias is lower than 0.2, round and then derived the median estimates of AP across 126 or (3) the number of validation checks failure reached 100 times LOO rounds to define the AP maps and their uncertainties across (i.e., the predicting bias for validation sets doesn’t decrease for the European continent. 100 times). Criterion (1) and (2) helps to avoid the over-fitting issue. Criterion (3) helps to avoid additional simulations without Attribution of AP Pattern to Drivers Across any improvement of BPN model. Europe ML models like BPN are able to extract hidden, complex, and Regression Trees generally subtle, non-linear relationships from scattered data BPN models cannot rate the importance of predictors for and multiple predictors and use these relationships for up- accounting for differences in observed AP, so we also used RTs scaling. The up-scaled AP maps are necessary but difficult to be to build regression models. The input predictors, training and interpreted by ML models due to the numbers of predictors and testing subset selection, leave-one-out validation and best model their interactions. We thus applied partial correlation analysis, a selection for the RT were the same as those used in the BPN simple diagnostic analysis, on modeled maps of AP to gain a first method. The response data at each binary split were grouped into order understanding of the most important predictors emerging two descendant nodes to maximize homogeneity, and the best from the model results in different regions. binary split was selected. The trees were grown to their maximum The spatial partial correlation analysis for Europe was sizes following the tree template. We used an RT tree template conducted at spatial resolution of half degree (∼50 km). To do where the minimal number of parents was 6, and the minimal so, predictors were firstly resampled into half degree by using the number of observations at the terminal nodes of the trees was 3. same methods as described above. Then we performed the partial correlation analysis between AP and predictors (only numerical predictors) using a spatial moving window of 4.5 × 4.5◦ . This Leave-One-Out Validation Method and analysis partly avoided the collinearity among the predictors but still failed in separating effects of strongly collinear factors BPN and RT Model Selection (e.g., NPP and climatic variables). Therefore, we performed (1) The leave-one-out (LOO) cross validation method is best suited a partial correlation analysis between AP and all abiotic factors, for model calibration and evaluation when observational data and (2) a partial correlation analysis between AP, soil conditions are scarce (Allen, 1974). LOO can test the reproducibility of and NPP without climatic variables. The median estimate of AP prediction models using independent data (Cawley, 2006). Out of and the predictor metrics for Europe were first resampled to the 126 European sites (Data set A), 125 observations were used 0.5◦ . We thus selected 81 pixels surrounding each 0.5 grid cell to train and select models and the remaining observation to verify for calculating the coefficients of partial correlation for AP and the model. This procedure was repeated 126 times by excluding the predictors. each observation once. Within the BPN and RT models training and selection, the 126 observations were randomly separated into training (106 RESULTS sites, 85%) and testing (19 sites, 15%) data sets (Figure 1). The BPN models were trained 5,000 times using different training Model Performance Across European Sites subsets, and the RT models were trained 500 times. We used and Along Spatial Gradients R2 and RMSE to assess the performance of the BPN and RT We assessed the performance of the BPN and RT models for each models on the training and testing data sets for each train-test of the 126 European sites (Data set A) by LOO cross-validation. possibility. We removed BPN and RT models with training R2 Excluding nine outliers with absolute bias >20 µmol g−1 h−1 <0.6, testing R2 <0.5, training RMSE >15 µmol g−1 h−1 and and relative bias >50% (Supporting Information Appendix S4, testing RMSE >10 µmol g−1 h−1 . Then, we selected the 25 Figures S4, S5), the BPN models explained 58% of AP variability models with the highest test-R2 for observed AP < 25 µmol g−1 across all sites, with an RMSE of 6.83 µmol g−1 h−1 (Figure 3, h−1 and lowest RMSE for observed AP >25 µmol g−1 h−1 . This Table S2). The BPN models reproduced 93% of the variance selection algorithm ensured that the BNP and RT can accurately of AP for the remaining 117 sites. The observations of AP for reproduce AP in both the training and testing data sets. 68 sites were within the uncertainty quantile range of 10–90% Frontiers in Big Data | www.frontiersin.org 4 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 2 | Spatial pattern of the 12 predictors: (A) soil pH (soilpH), (B) soil organic carbon content (OC), (C) soil total nitrogen content (TN), (D) soil total phosphorus content (TP), (E) clay content (Clay), (F) soil labile organic phosphorus content (Po), (G) net primary productivity (NPP), (H) mean annual temperature (MAT), (I) amplitude of yearly temperature (AMP), (J) mean annual precipitation (MAP), (K) vegetation cover (VegType), and (L) soil type (SoilType). of AP simulated with the ensemble of 25 BPN models, and France, Germany, Poland and Sweden. Median estimates of 32 site observations were within the 25–75% quantile range modeled AP along this transect gradually decrease from 8.5 µmol (Figure 3). The BPN model generally underestimated AP, but g−1 h−1 in Spain to 5.6 µmol g−1 h−1 for Germany, which is overestimated AP at sites with low activity (<20 µmol g−1 roughly consistent with the observations (decrease from 11.8 to h−1 ), while AP at sites with intermediate and high observed 2.3 µmol g−1 h−1 ). Modeled AP for sites in Poland and Sweden activities (>20 µmol g−1 h−1 ) tended to be systematically were all higher than observed values, though. Transect 2 (across underestimated. The BPN models outperformed the RT models England and France to Italy) shows observed values going from for LOO validation (Table S2). With excluding 13 outliers, RT 7.8 µmol g−1 h−1 (6.6–9.9 µmol g−1 h−1 ) in England to 4.2 µmol models explained 48.7% of the spatial variance of AP over all g−1 h−1 (0.4–9.8 µmol g−1 h−1 ) in Italy. Predicted AP along measurements, accompanied by a higher RMSE of 7.85 µmol g−1 Transect 2 give a smaller decrease than observed, going from 7.4 h−1 (Figure S6, Table S2). µmol g−1 h−1 (6.1–9.3 µmol g−1 h−1 ) to 7.2 µmol g−1 h−1 (5.5– Three transects (Data set B) were selected to determine if 8.8 µmol g−1 h−1 ). Modeled AP along Transect 3 (from Western the BPN model could reproduce large scale observed spatial and Eastern Spain, Italy and Turkey) decreased from 8.5 µmol gradients (Figures 4C–E). Transect 1 spanned across Spain, g−1 h−1 (5.1–10.8 µmol g−1 h−1 ) for Western Spain to 5.4 µmol Frontiers in Big Data | www.frontiersin.org 5 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 3 | Performance of the BPN models on 126 European sites following the leave-one-out framework. (A) Predicted vs. observed acid phosphatase activity (AP) for all sites. (B) Predicted vs. observed AP with excluding seven outliers (gray circles). The dashed lines indicate the 1:1 lines. The blue lines indicate the regression lines between predicted and observed AP. The blue and gray error bars indicate the 25 and 75% quantiles of predicted AP, and the orange error bars indicate the 10 and 90% quantiles of predicted AP. The nine outliers are marked in gray: Site A and D are located in central England, Site B is in Northern Italy, Site C in Eastern Spain, Site E in Poland, Site F in Norway, and Site G in Western Spain. g−1 h−1 for Turkey, a smaller decrease than the observations. variables (Figures 6G–I,K) which reflects strong collinearity This bias could be mainly attributed to the underestimation of between NPP and climatic variables (Figure S9). modeled AP for sites with high values in Spain. Modeled AP is found to be negatively correlated with soil In summary, The BPN models are found in reproduce the pH at high soil pH values (6–9), which is common for Europe, spatial variance of AP for Europe and the direction of the AP especially in the Southern and the Western of the domain gradients across ecological zones and climate, despite negative (Figure 6A). Soil nutrient contents also partly accounted for the or positive biases compared to individual site-scale data points. spatial patterns of modeled AP. AP is found to be positively Median BPN modeled estimate of AP across Europe increases correlated with soil organic C content (OC) (Figures 6B–D,F). from north to south, and also from east to west (Figure 4A). More specifically, AP and OC are significantly positively AP was highest for Northern England and Western Spain (10–40 correlated (p < 0.1) in Southern and Western Europe but µmol g−1 h−1 ), and lowest for Northern and Eastern Europe (<2 not significantly negatively correlated for Northern and central µmol g−1 h−1 ) (Figure 4). The uncertainty in the AP estimates Europe. AP and TN are strongly positively correlated with AP (the difference between 75 and 25% quantiles; Figure S7) was low in Southern and the Western Europe. AP is significantly (p for most of Europe (<3 µmol g−1 h−1 ) except for England and < 0.1) positively correlated with TP but negatively correlated Northern Europe (>3 µmol g−1 h−1 ). with Po in aridity area (P/PET < 0.7) (Figure 6; Figure S11), while AP is significantly (p < 0.1) negatively correlated with TP Relationship Between Modeled AP and but positively correlated with Po in wetter region (P/PET > 1; Climate and Soil Conditions Figure 6; Figure S11). The partial correlation analysis indicated that correlations Modeled AP is found to be higher in highly weathered between gridded TP and Po and gridded climatic factors (MAT, soils (Ultisols and Oxisols, 7.8–8.3 µmol g−1 h−1 ) compared AMP, and MAP) with gridded BPN estimates of AP for Europe to slightly and moderately weathered soils (1.2–7.2 µmol g−1 were stronger than with the other factors (Figure 5). The h−1 ) (Figure S8A). AP was the highest (8.9 µmol g−1 h−1 ) proportions of land areas with significant (p < 0.1) correlations for Evergreen Broadleaf Forests and the lowest for Evergreen between AP and MAT, AMP and MAP were 76, 79, and 71%, Needleleaf Forests among all of the biomes (Figure S8B). respectively. AP was generally correlated with both MAT and MAP (Figures 6G–I), being positively associated with spatially Upscaling of Site Data to a Global Map of increasing MAT, and negatively correlated with AMP. AP When excluding climatic variables from partial correlation We used the BPN models trained on European data in analysis, NPP has a significant correlation (p < 0.1) with AP a first attempt to scale up AP to the global scale. The for 58% of land areas. Significant positive correlations between modeled global spatial pattern of AP show a clear latitudinal NPP and AP are mainly found in Southern and Western Europe, gradient of increasing AP from high latitudes to the equator while significant negative correlations are in Sweden and Finland (Figure 7). Modeled AP is also much higher for coastal (Figure 6K). The emerging relationship between modeled AP than inland areas. The global up-scaled distribution of AP and NPP is consistent with those found between AP and climatic is however more uncertain than the European distribution. Frontiers in Big Data | www.frontiersin.org 6 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 4 | Spatial pattern of median estimates of acid phosphatase activity (AP) from the BPN models for Europe (A). The observed AP are also shown. (B) Uncertainties of AP defined as the difference between the 75 and 25% quantiles of predicted AP from the BPN models. (C–E) show the changes of observed and predicted AP along three transects: (C) Spain (SP), France (FR), Germany (GE), Poland (PO), and Sweden (SW); (D) England (EN), France (FR), and Italy (IT); and (E) Western Spain (west of 5◦ W; SP1), Eastern Spain (east of 5◦ W; SP2), Italy (IT), and Turkey (TR). FIGURE 5 | Proportions of vegetated land in Europe controlled by (A) climatic (red) and soil correlation analysis between predicted acid phosphatase activity (AP) and all numerical abiotic factors. (B) Was based on the partial correlation analysis between AP and NPP and soil factors. Pixels with p > 0.1 are excluded. The soil factors are: soil pH (SoilpH), soil organic carbon content (OC), soil total nitrogen content (TN), soil total phosphorus content (TP), soil clay content (Clay), soil labile organic phosphorus content (Po). The climatic factors are: mean annual temperature (MAT), amplitude of yearly temperature (AMP), and mean annual precipitation (MAP). For example, AP for about 48% of the vegetated land is controlled by soil pH, with partial correlation coefficients between AP and soil pH > 0.2. The BPN models explained only a small fraction (<15%) of set C) and tropical regions (Data set D) (Figures S12A–D). The AP (Supporting Information Appendix S5) variance over all low prediction capacity for these sites may be partly attributed measurements for both temperate sites outsides Europe (Data to the low accuracy of extracted gridded TN from ISRIC-WISE Frontiers in Big Data | www.frontiersin.org 7 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 6 | Partial correlations between the median estimates of acid phosphatase activity (AP) and 9 abiotic factors: (A) soil pH, (B) soil organic carbon, (C) soil total nitrogen content, (D) soil total phosphorus content, (E) soil clay content, (F) soil labile organic phosphorus content, (G) mean annual temperature, (H) amplitude of yearly temperature, and (I) mean annual precipitation. (J,K) Are the partial correlations between AP and NPP when including climatic variables and excluding climatic variables. Pixels with p > 0.1 were excluded. Frequency histograms showing the areal coverages correspond to positive and negative correlations, estimated as proportions of the total study area, are provided for each map (insets). (Supporting Information Appendix S1). Excluding the sites continent where the training data originate from, and still has with gap-filled TN, the BPN models explained 34% of AP some predictive ability at cross validation sites that originate from variability for the rest of temperate sites outside Europe other temperate regions, mainly USA and China (Figure 7). (Figures S12B,C). Besides, MAT of many sites fall out of the range of the Europe training sites, which also lead to the poor performance of our BNP models on sites outside Europe. DISCUSSION Excluding sites with gap-filled TN and sites with MAT falling out the ranges between 5 and 20◦ C, the BPN models explained 42% Extrapolation of AP From Site to Regional of AP variability for temperate sites outside Europe within the and Global Scales MAT range of the training sites (5◦ C < MAT < 20◦ C) (Figure 3). This study is a first attempt to extrapolate AP site measurements This indicates that the BPN can predict spatial gradients of from Europe to regional and global scales using ML methods. AP with a reasonable accuracy (R2 = 0.58) over the European The BPN method successfully reproduced the spatial variance Frontiers in Big Data | www.frontiersin.org 8 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 7 | Global pattern of AP from the BPN models (median). “x” indicated area with mean annual temperature (MAT) < 5◦ C or MAT > 20◦ C which falling out of the ranges of BPN built based on Europe sites. The observed AP are shown as dots. of AP across Europe. The performance and reproducibility of from manipulation experiments (soil warming and irrigation). BPN models strongly depend on the quantity, quality, and Warming had a positive effect on AP in experimental studies for representativeness of the training data set. We attribute the both temperate grassland soil (0–10 cm depth; Zhou et al., 2013) failure in extrapolation to global scale (Figure 8) to an insufficient and Mediterranean shrubland (Sardans et al., 2006). A possible number of representative training AP measurements, especially mechanism is that increasing temperature can accelerate plant for tropical soils and climate. growth and P acquisition by plants, thus reducing soil P content AP activity in tropical regions varies widely from 0 and stimulating AP production (Sardans et al., 2006). The to 80 µmol g−1 h−1 . The BPN model trained across depletion of soil P can stimulate the production of AP by both Europe was not able to capture this large range. A low plants and soil microbes (Clarholm, 1993; Olander and Vitousek, representativeness of observations in Europe for vegetation types 2000). Drought experiments in the Mediterranean reported no and climatic conditions outside Europe can partly account for response (Sardans and Peñuelas, 2005) or reduction in AP this low prediction accuracy. The more homogenous climatic activity (Sardans et al., 2006) under conditions of low soil-water conditions in the tropics than the extra-tropics may increase the content, associated with low water availability and plant growth importance of soil conditions in accounting for AP variation and ultimately with a reduction in the P demand of plants and (Supporting Information Appendix S6). Soil properties in the microbes. Zhou et al. (2013) reported that increased precipitation tropics are very heterogeneous and were shown to drive large significantly increased soil AP in a temperate grassland at depths scale patterns of ecosystem properties (Quesada et al., 2012). of 10–20 cm. Measurements based on spatial sampling have also More AP measurements in tropical regions would be needed to indicated that a higher availability of soil water leads to a larger help us understand the differences between soil types and P-use microbial biomass (Brockett et al., 2012), and hence increased strategies of ecosystems in this region. activities of soil extracellular enzymes. Increasing precipitation should also promote substrate diffusion, which brings more Drivers of the Spatial Patterns of AP in soil NO3 -N and NH4 -N to the subsurface and increases AP Europe (Zhou et al., 2013). The assessment of the major factors influencing potential AP help us to understand (1) the controls on AP across Europe and (2) Relationships Between Modeled AP and Soil the failure of our model to extrapolate to a global scale. Simple Variables regression analyses were applied in previous studies to investigate Nutrient availability, pH, organic C, and clay content are also the factors influencing AP (Sinsabaugh et al., 2008; Huang et al., found to play important roles in the modeled spatial variance of 2013; Margalef et al., 2017). Such approaches, however, cannot AP from west to east Europe (Figures 6A–F). Modeled AP was manage collinearity among predictors nor account for non-linear positively correlated with SOC. The substrate for AP, soil organic responses to different variables (Ma et al., 2017). We applied P, is tightly linked to SOC (Tipping et al., 2016), but not included partial correlation diagnostic to the output of the BPN model, in our analysis due to lack of data. Therefore, the correlation which partly overcomes these shortcomings in separating the between AP and SOC can be interpreted as a control of substrate effect of each predictor on AP for different climatic zones and on AP. This finding was consistent with results from a meta- states of soil nutrients. analysis of soil phosphatase activity and transect measurements (Sinsabaugh et al., 2008; Huang et al., 2013). Relationships Between Modeled AP and Climatic TN was also positively correlated with modeled AP for most Variables of Europe (Figure 6C). This finding is consistent with results Modeled AP is found to be positively correlated with MAT and by Margalef et al. (2017) and a meta-analysis of N-fertilization MAP for most of Europe (Figures 6G,I) in line with results experiments by Marklein and Houlton (2012). Higher soil N Frontiers in Big Data | www.frontiersin.org 9 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity FIGURE 8 | Performance of the back-propagation networks on temperate sites (excluding Europe) (A) and tropical sites (B). Temperate and tropical sites with original TN and 5◦ C < MAT < 20◦ C were reproduced using BPNs and were then compared to observations. The solid lines indicate the regression lines between predicted vs. observed AP. The dashed lines indicate the 1:1 lines. The error bars indicate the 10 and 90% quantiles of predicted AP. availability is supposed to strongly increase AP activity in both from most experimental studies, which provided evidence of a plant roots and bulk soils across a wide range of ecosystems positive association between AP and environmental resources (Marklein and Houlton, 2012). and conditions (water, temperature, and nutrients). Our analysis In contrast to the monotonic behavior of AP with TN and showed that MODIS-NPP, which covaries with climate, is also climatic variables, soil TP and Po were non-monotonically statistically associated with modeled AP (Figure 6K). Generally, correlated with modeled AP for Europe (Figures 6D,F). AP in our results are consistent with the mechanism that adequate wet regions (P/PET > 0.7) was negatively correlated with TP but resources (temperature, water and nutrients) can sustain a positively correlated with Po. This could indicate that biological higher production and P demand for plants and microbes, thus mineralization of Po via AP is an essential strategy to maintain larger AP activity to increase mineralization rates (Margalef the bioavailable P pool under humid conditions (Izquierdo et al., et al., 2017). However, the co-relationship between AP and 2013). AP in aridity area (P/PET < 0.7) is positively correlated NPP is weaker than that with climatic variables and TN with TP but negatively correlated with Po (Figure 6). Compared (Figures 6G–I,K), which indicating that the AP pattern is not with humid soils, arid soils typically contain lower SOM, thus only determined by P demand by biota but also by climate and the importance of minerals as source of P for biota is suggested soil nutrients. to become more important (Lajtha and Schlesinger, 1988; Cross and Schlesinger, 1995; Delgado-Baquerizo et al., 2013). Our Relationships Between Modeled AP and Vegetation results show a lower investment into AP under semi-dry to and Soil Types dry conditions providing support for such shifts in nutritional Modeled AP differed among the soil types, with higher activities strategies on large spatial scales (e.g., Feng et al., 2016). in Oxisols and Ultisols (Figure S8) in line with earlier findings Soil pH was negatively correlated with modeled AP across (Acosta-Martinez et al., 2007). These soils are commonly poor in Europe (Figure 6A), in agreement with the AP measurements P and thus a high recycling efficiency of organic P is an adaption based on acidic soils (Dick et al., 2000). Soils with high pHs to these conditions (Vance et al., 2003). We found that AP are generally associated with fewer amino acid functional groups was higher in evergreen broadleaf forest than deciduous forests that are essential for catalysis (Dick et al., 2000). Clay content (Figure S8), consistent with the finding by Redel et al. (2008) had a non-linear but weak relationship with modeled AP across of higher litter accumulation and lower P cycling in deciduous Europe: increasing clay content increased AP in soils with low forests. Margalef et al. (2017) found that soil phosphatase was clay content (<10%) but decreased AP in soils with high clay significantly higher in angiosperm than gymnosperm forests, content (>10%) (Figure 6E). Clay particles are characterized likely because gymnosperm tend to grow on soils which by many reactive binding sites for P and AP leading to the have a lower biological activity as an adaptation to harsh deactivation of the latter (Rao et al., 2000). Thus, other strategies conditions and thus have the associated disadvantage of low than AP, like extraction of chelating agents and acids, might be N-use efficiency. advantageous in clay rich soils. Implication and Uncertainties Relationships Between Modeled AP and Net Primary The BPN models successfully extrapolated the site AP Production measurements on the Europe region. However, there are The relationships between abiotic (climate, soils) factors and some inevitable uncertainties due to the gap-filling of input AP generally are in qualitative agreement with the results predictors which were not reported in the original literature. Frontiers in Big Data | www.frontiersin.org 10 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity The extracted values from gridded datasets are biased to some global models (e.g., Goll et al., 2012). Further, we identified extent, e.g., due to the coarse spatial resolution of predictor the need for more measurements in tropical, arid, boreal, or maps as some variables like TN, TP, and Clay vary strongly alpine ecosystems, so as to be able to extend our approach to within very short (cm-m) distance (Lechowicz and Bell, 1991; global scale. Iqbal et al., 2005). As a consequence, the covariation between observed and extracted values of TN, TP, and Clay are not very DATA AVAILABILITY STATEMENT strong (Figure S1). However, non-linear BPN models are considered to be The datasets of predictors used in this study are publicly more sensitive to the orders of predictor values rather than available and can be accessed by the links we provided linear variation among predictor values. This advantage of BPN in Supporting Information. The code for BPN and RT compared to linear models reduces the influence of biases in can be accessed on Github: https://github.com/yansun8/AP_ gap-filling on AP predictions. The spearman rank correlations extropolation. Other data are available upon request from between observed and extracted values for Europe sites are found the authors. to be significant for all variables except for Clay (Figure S1). Thus, the data quality seems in general sufficient for the BPN model. Still, the BPN model might underestimate the driving AUTHOR CONTRIBUTIONS effect of Clay on AP pattern due to the large biases. High-quality YS, DG, and SP designed the research. YS and DG wrote the gridded datasets especially for soil clay content are required, to paper. YS did the extrapolation of AP and analyzed the data. OM reduce the uncertainty in the extrapolations based on machine and DA collected the dataset of site measurements. PC, JS, and JP learning methods. helped revise this paper. CONCLUSIONS FUNDING This study is the first quantification of AP on continental scale by extrapolating the site observations of AP to a This work was supported by the European Research Council spatial explicit map of AP which can be used to benchmark Synergy project SyG-2013-610028 IMBALANCE-P. land surface models for temperate regions. This study also provided the comprehensive overview of the factors driving SUPPLEMENTARY MATERIAL the spatial variations of phosphatase activity on continental scale. The extent to which plants can upregulate P recycling The Supplementary Material for this article can be found online under future climate changes could be limited by climatic at: https://www.frontiersin.org/articles/10.3389/fdata.2019. and soil factors which are currently not accounted for in 00051/full#supplementary-material REFERENCES Clarholm, M. (1993). Microbial biomass P, labile P, and acid phosphatase activity in the humus layer of a spruce forest, after repeated additions of fertilizers. Biol. Acosta-Martinez, V., Cruz, L., Sotomayor-Ramirez, D., and Pérez-Alegría, L. Fertil. Soils 16, 287–292. doi: 10.1007/BF00369306 (2007). Enzyme activities as affected by soil properties and land use in a tropical Cross, A. F., and Schlesinger, W. H. (1995). A literature review and watershed. Appl. Soil Ecol. 35, 35–45. doi: 10.1016/j.apsoil.2006.05.012 evaluation of the Hedley fractionation: applications to the biogeochemical Allen, L. H. (1974). Model of light penetration into a wide-row crop 1. Agron. J. 66, cycle of soil phosphorus in natural ecosystems. Geoderma 64, 197–214. 41–47. doi: 10.2134/agronj1974.00021962006600010011x doi: 10.1016/0016-7061(94)00023-4 Ballabio, C., Lugato, E., Fernández-Ugalde, O., Orgiazzi, A., Jones, A., Delgado-Baquerizo, M., Maestre, F. T., Gallardo, A., Bowker, M. A., Wallenstein, Borrelli, P., et al. (2019). Mapping LUCAS topsoil chemical properties at M. D., Quero, J. L., et al. (2013). Decoupling of soil nutrient cycles as a European scale using Gaussian process regression. Geoderma 355:113912. function of aridity in global drylands. Nature 502, 672–676. doi: 10.1038/nature doi: 10.1016/j.geoderma.2019.113912 12670 Ballabio, C., Panagos, P., and Monatanarella, L. (2016). Mapping topsoil physical Dick, W. A., Cheng, L., and Wang, P. (2000). Soil acid and alkaline phosphatase properties at European scale using the LUCAS database. Geoderma 261, activity as pH adjustment indicators. Soil Biol. Biochem. 32, 1915–1919. 110–123. doi: 10.1016/j.geoderma.2015.07.006 doi: 10.1016/S0038-0717(00)00166-8 Brockett, B. F. T., Prescott, C. E., and Grayston, S. J. (2012). Soil moisture is the Eivazi, F., and Tabatabai, M. A. (1977). Phosphatases in soils. Soil Biol. Biochem. 9, major factor influencing microbial community structure and enzyme activities 167–172. doi: 10.1016/0038-0717(77)90070-0 across seven biogeoclimatic zones in western Canada. Soil Biol. Biochem. 44, Ellsworth, D. S., Anderson, I. C., Crous, K. Y., Cooke, J., Drake, J. E., 9–20. doi: 10.1016/j.soilbio.2011.09.003 Gherlenda, A. N., et al. (2017). Elevated CO2 does not increase eucalypt forest Caldwell, B. A. (2005). Enzyme activities as a component of soil biodiversity: a productivity on a low-phosphorus soil. Nat. Clim. Chang. 7:279. doi: 10.1038/ review. Pedobiologia 49, 637–644. doi: 10.1016/j.pedobi.2005.06.003 nclimate3235 Cawley, G. C. (2006). “Leave-one-out cross-validation based model selection Feng, J., Turner, B. L., Lü, X., Chen, Z., Wei, K., Tian, J., et al. (2016). Phosphorus criteria for weighted LS-SVMs,” in The 2006 IEEE International Joint transformations along a large-scale climosequence in arid and semiarid Conference on Neural Network Proceedings (Vancouver, BC), 1661–1668. grasslands of northern China. Glob. Biogeochem. Cycles 30, 1264–1275. doi: 10.1109/IJCNN.2006.246634 doi: 10.1002/2015GB005331 Frontiers in Big Data | www.frontiersin.org 11 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity Fick, S. E., and Hijmans, R. J. (2017). WorldClim 2: new 1-km spatial resolution Marklein, A. R., and Houlton, B. Z. (2012). Nitrogen inputs accelerate phosphorus climate surfaces for global land areas. Int. J. Climatol. 37, 4302–4315. cycling rates across a wide variety of terrestrial ecosystems. New Phytol. 193, doi: 10.1002/joc.5086 696–704. doi: 10.1111/j.1469-8137.2011.03967.x Friedl, M. A., Sulla-Menashe, D., Tan, B., Schneider, A., Ramankutty, N., Sibley, Mason, C., Twomey, J., Wright, D., and Whitman, L. (2017). Predicting A., et al. (2010). MODIS Collection 5 global land cover: algorithm refinements engineering student attrition risk using a probabilistic neural network and and characterization of new datasets. Rem. Sens. Environ. 114, 168–182. comparing results with a backpropagation neural network and logistic doi: 10.1016/j.rse.2009.08.016 regression. Res. High Educ. 59, 382–400. doi: 10.1007/s11162-017-9473-z Gill, A. L., and Finzi, A. C. (2016). Belowground carbon flux links biogeochemical Nannipieri, P., Giagnoni, L., Landi, L., and Renella, G. (2011). “Role of cycles and resource-use efficiency at the global scale. Ecol. Lett. 19, 1419–1428. phosphatase enzymes in soil,” in Phosphorus in Action, eds E. Bünemann, doi: 10.1111/ele.12690 A. Oberson, and E. Frossard (Berlin; Heidelberg: Springer), 215–243. Goll, D. S., Brovkin, V., Parida, B. R., Reick, C. H., Kattge, J., Reich, P. B., et al. doi: 10.1007/978-3-642-15271-9_9 (2012). Nutrient limitation reduces land carbon uptake in simulations with a Olander, L. P., and Vitousek, P. M. (2000). Regulation of soil phosphatase model of combined carbon, nitrogen and phosphorus cycling. Biogeosciences 9, and chitinase activity by N and P availability. Biogeochemistry 49, 175–191. 3547–3569. doi: 10.5194/bg-9-3547-2012 doi: 10.1023/A:1006316117817 Golubiewski, N. E. (2006). Urbanization increases grassland carbon pools: Papale, D., and Valentini, R. (2003). A new assessment of European forests carbon effects of landscaping in Colorado’s front range. Ecol. Appl. 16, 555–571. exchanges by eddy fluxes and artificial neural network spatialization. Glob. doi: 10.1890/1051-0761(2006)016[0555:uigcpe]2.0.co;2 Change Biol. 9, 525–535. doi: 10.1046/j.1365-2486.2003.00609.x Hengl, T., Mendes de Jesus, J., Heuvelink, G. B., Ruiperez Gonzalez, M., Quesada, C. A., Phillips, O. L., Schwarz, M., Czimczik, C. I., Baker, T. R., Kilibarda, M., Blagotić, A., et al. (2017). SoilGrids250m: Global gridded Patiño, S., et al. (2012). Basin-wide variations in Amazon forest structure and soil information based on machine learning. PLoS ONE 12:e0169748. function are mediated by both soils and climate. Biogeosciences 9, 2203–2246. doi: 10.1371/journal.pone.0169748 doi: 10.5194/bg-9-2203-2012 Hogan, E. J., Minnullina, G., Smith, R. I., and Crittenden, P. D. (2010). Rao, M. A., Violante, A., and Gianfreda, L. (2000). Interaction of Effects of nitrogen enrichment on phosphatase activity and nitrogen: acid phosphatase with clays, organic molecules and organo-mineral phosphorus relationships in Cladonia portentosa. New Phytol. 186, 911–925. complexes: Kinetics and stability. Soil Biol. Biochem. 32, 1007–1014. doi: 10.1111/j.1469-8137.2010.03222.x doi: 10.1016/S0038-0717(00)00010-9 Huang, W., Liu, J., Wang, Y. P., Zhou, G., Han, T., and Li, Y. (2013). Increasing Redel, Y., Rubio, R., Godoy, R., and Borie, F. (2008). Phosphorus fractions and phosphorus limitation along three successional forests in southern China. Plant phosphatase activity in an Andisol under different forest ecosystems. Geoderma Soil 364, 181–191. doi: 10.1007/s11104-012-1355-8 145, 216–221. doi: 10.1016/j.geoderma.2008.03.007 Huang, W., Liu, J., Zhou, G., Zhang, D., and Deng, Q. (2011). Effects of Reed, S. C., Yang, X., and Thornton, P. E. (2015). Incorporating phosphorus cycling precipitation on soil acid phosphatase activity in three successional forests in into global modeling efforts: a worthwhile, tractable endeavor. New Phytol. 208, southern China. Biogeosciences 8, 1901–1910. doi: 10.5194/bg-8-1901-2011 324–329. doi: 10.1111/nph.13521 Iqbal, J., Thomasson, J. A., Jenkins, J. N., Owens, P. R., and Whisler, F. D. (2005). Riedmiller, M., and Braun, H. (1993). “A direct adaptive method for faster Spatial variability analysis of soil physical properties of alluvial soils. Soil Sci. backpropagation learning: the RPROP algorithm,” in Proceedings of the IJCNN Soc. Am. J. 69, 1338–1350. doi: 10.2136/sssaj2004.0154 (San Francisco, CA: IEEE Press), 586–591. Izquierdo, J. E., Houlton, B. Z., and Huysen, T. L. V. (2013). Evidence for Running, S. W., Nemani, R. R., Heinsch, F. A., Zhao, M., Reeves, M., progressive phosphorus limitation over long-term ecosystem development: and Hashimoto, H. (2004). A continuous satellite-derived measure examination of a biogeochemical paradigm. Plant Soil 367, 135–147. of global terrestrial primary production. Bioscience 54, 547–560. doi: 10.1007/s11104-013-1683-3 doi: 10.1641/0006-3568(2004)054[0547:ACSMOG]2.0.CO;2 Keenan, T. F., Davidson, E., Moffat, A. M., Munger, W., and Richardson, Sardans, J., and Peñuelas, J. (2005). Drought decreases soil enzyme A. D. (2012). Using model-data fusion to interpret past trends, and activity in a Mediterranean. Soil Biol. Biochem. 37, 455–461. quantify uncertainties in future projections, of terrestrial ecosystem carbon doi: 10.1016/j.soilbio.2004.08.004 cycling. Glob. Change Biol. 18, 2555–2569. doi: 10.1111/j.1365-2486.2012. Sardans, J., Peñuelas, J., and Estiarte, M. (2006). Warming and drought alter soil 02684.x phosphatase activity and soil P availability in a Mediterranean shrubland. Plant Kelley, H. J. (1960). Gradient theory of optimal flight paths. ARS J. 30, 947–954. Soil 289, 227–238. doi: 10.1007/s11104-006-9131-2 doi: 10.2514/8.5282 Shangguan, W., Dai, Y., Duan, Q., Liu, B., and Yuan, H. (2014). A global soil Kitayama, K. (2013). The activities of soil and root acid phosphatase in the nine data set for earth system modeling. J. Adv. Model. Earth. Syst. 6, 249–263. tropical rain forests that differ in phosphorus availability on Mount Kinabalu, doi: 10.1002/2013MS000293 Borneo. Plant Soil 367, 215–224. doi: 10.1007/s11104-013-1624-1 Sinsabaugh, R. L., Lauber, C. L., Weintraub, M. N., Ahmed, B., Allison, S. D., Krämer, S., and Green, D. M. (2000). Acid and alkaline phosphatase dynamics Crenshaw, C., et al. (2008). Stoichiometry of soil enzyme activity at global scale. and their relationship to soil microclimate in a semiarid woodland. Soil Biol. Ecol. Lett. 11, 1252–1264. doi: 10.1111/j.1461-0248.2008.01245.x Biochem. 32, 179–188. doi: 10.1016/S0038-0717(99)00140-6 Sun, Y., Peng, S., Goll, D. S., Ciais, P., Guenet, B., Guimberteau, M., et al. (2017). Lajtha, K., and Schlesinger, W. H. (1988). The biogeochemistry of phosphorus Diagnosing phosphorus limitations in natural terrestrial ecosystems in carbon cycling and phosphorus availability along a desert soil chronosequence. Ecology cycle models. Earth’s Fut. 5, 730–749. doi: 10.1002/2016EF000472 69, 24–39. doi: 10.2307/1943157 Tipping, E., Somerville, C. J., and Luster, J. (2016). The C: N: P: S Lechowicz, M. J., and Bell, G. (1991). The ecology and genetics of fitness in forest stoichiometry of soil organic matter. Biogeochemistry 130, 117–131. plants. II. Microspatial heterogeneity of the edaphic environment. J. Ecol. 79, doi: 10.1007/s10533-016-0247-z 687–696. doi: 10.2307/2260661 Turner, D. P., Ritts, W. D., Cohen, W. B., Gower, S. T., Running, S. W., Zhao, Ma, Z., Liu, H., Mi, Z., Zhang, Z., Wang, Y., Xu, W., et al. (2017). Climate warming M., et al. (2006). Evaluation of MODIS NPP and GPP products across multiple reduces the temporal stability of plant community biomass production. Nat. biomes. Rem. Sens. Environ. 102, 282–292. doi: 10.1016/j.rse.2006.02.017 Commun. 8, 1–7. doi: 10.1038/ncomms15378 Vance, C. P., Uhde-Stone, C., and Allan, D. L. (2003). Phosphorus acquisition and Maistry, P. M., Muasya, A. M., Valentine, A. J., and Chimphango, S. B. M. (2015). use: critical adaptations by plants for securing a nonrenewable resource. New Increasing nitrogen supply stimulates phosphorus acquisition mechanisms Phytol. 157, 423–447. doi: 10.1046/j.1469-8137.2003.00695.x in the fynbos species Aspalathus linearis. Funct. Plant Biol. 42, 52–62. Vitousek, P. M., Porder, S., Houlton, B. Z., and Chadwick, O. (2010). Terrestrial doi: 10.1071/FP14100 phosphorus limitation: mechanisms, implications, and nitrogen–phosphorus Margalef, O., Sardans, J., Fernández-Martínez, M., Molowny-Horas, R., Janssens, interactions. Ecol. Appl. 20, 5–15. doi: 10.1890/08-0127.1 I. A., Ciais, P., et al. (2017). Global patterns of phosphatase activity in natural Walker, T., and Syers, J. (1976). The fate of phosphorus during pedogenesis. soils. Sci. Rep. 7, 1–13. doi: 10.1038/s41598-017-01418-8 Geoderma 15, 1–19. doi: 10.1016/0016-7061(76)90066-5 Frontiers in Big Data | www.frontiersin.org 12 January 2020 | Volume 2 | Article 51 Sun et al. Spatial Pattern of Acid Phosphatase Activity Wang, R., Balkanski, Y., Boucher, O., Ciais, P., Peñuelas, J., and Tao, S. (2015). Zhao, M., Heinsch, F. A., Nemani, R. R., and Running, S. W. (2005). Significant contribution of combustion-related emissions to the atmospheric Improvements of the MODIS terrestrial gross and net primary production phosphorus budget. Nat. Geosci. 8:48. doi: 10.1038/ngeo2324 global data set. Rem. Sens. Environ. 95, 164–176. doi: 10.1016/j.rse.2004. Were, K., Bui, D. T., Dick, Ø. B., and Singh, B. R. (2015). A comparative assessment 12.011 of support vector regression, artificial neural networks, and random forests Zhou, X., Chen, C., Wang, Y., Xu, Z., Han, H., Li, L., et al. (2013). Warming for predicting and mapping soil organic carbon stocks across an Afromontane and increased precipitation have differential effects on soil extracellular landscape. Ecol. Indic. 52, 394–403. doi: 10.1016/j.ecolind.2014.12.028 enzyme activities in a temperate grassland. Sci. Total. Environ. 444, 552–558. White, P. J., and Hammond, J. P. (2008). “Phosphorus nutrition of doi: 10.1016/j.scitotenv.2012.12.023 terrestrial plants,” in The Ecophysiology of Plant-Phosphorus Interactions, eds P. J. White and J. P. Hammond (Dordrecht: Springer), 51–81. Conflict of Interest: The authors declare that the research was conducted in the doi: 10.1007/978-1-4020-8435-5_4 absence of any commercial or financial relationships that could be construed as a Wiesmeier, M., Barthold, F., Blank, B., and Kögel-Knabner, I. (2011). potential conflict of interest. Digital mapping of soil organic matter stocks using random forest modeling in a semi-arid steppe ecosystem. Plant Soil 340, 7–24, 2011. Copyright © 2020 Sun, Goll, Ciais, Peng, Margalef, Asensio, Sardans and doi: 10.1007/s11104-010-0425-z Peñuelas. This is an open-access article distributed under the terms of the Yang, X., Post, W. M., Thornton, P. E., and Jain, A. (2013). The distribution of soil Creative Commons Attribution License (CC BY). The use, distribution or phosphorus for global biogeochemical modeling. Biogeosciences 10, 2525–2537. reproduction in other forums is permitted, provided the original author(s) doi: 10.5194/bg-10-2525-2013 and the copyright owner(s) are credited and that the original publication in Yang, X., Thornton, P. E., Ricciuto, D. M., and Post, W. M. (2014). The role of this journal is cited, in accordance with accepted academic practice. No use, phosphorus dynamics in tropical forests–a modeling study using CLM-CNP. distribution or reproduction is permitted which does not comply with these Biogeosciences 11, 1667–1681. doi: 10.5194/bg-11-1667-2014 terms. Frontiers in Big Data | www.frontiersin.org 13 January 2020 | Volume 2 | Article 51