Zoning for edaphoclimatic aptitude of Bambusa vulgaris and Dendrocalamus giganteus in Brazil1

1. Received: Jun. 27, 2018. Accepted: Sep. 18, 2018. Published: Apr. 04, 2019. DOI: 10.1590/1983-40632019v4953713. 2. Universidade Federal de Lavras, Departamento de Ciências Florestais, Lavras, MG, Brasil. Email/ORCID: karinarodrigues_nina@hotmail.com/0000-0003-0963-0562, lucas.amaral@dcf.ufla.br/0000-0001-5219-9179, rlgrisi@dcf.ufla.br/0000-0001-7347-5523, sbotelho@dcf.ufla.br/0000-0003-4178-465X. 3. Universidade Estadual de Mato Grosso do Sul, Aquidauana, MS, Brasil. Email/ORCID: patpompeu@hotmail.com/ 0000-0002-3140-3457. The use of renewable materials and more sustainable processes are essential to minimize the negative impacts of population growth. Bamboo is a viable alternative, because it is a product with potential properties for wood replacement. However, to achieve a good yield, it is necessary to identify suitable sites for the forestry of the species. This study aimed to forecast favorable areas for the best establishment and development of Bambusa vulgaris and Dendrocalamus giganteus in the Brazilian territory, functioning as a zoning for edaphoclimatic aptitude that can support decisionmaking. It was performed using georeferenced points of natural occurrence for the species, collected from a virtual herbarium. The environmental variables were extracted from the WorldClim, MOD16 and SoilGrids databases. In order to select the set of variables to be used for modeling each species, a correlation test was performed. From this result, ten environmental variables were used for each species, and non-correlated variables were chosen. A principal component analyses was applied to select the variables that better explained variability. For the distribution map, the Maxent algorithm was used, based on the principle of maximum entropy, evaluated by the AUC and TSS metrics. Bambusa vulgaris presents a large area for adaptation in the Brazilian territory, unlike the Dendrocalamus giganteus species, which showed a limited aptitude area.

The use of renewable materials and more sustainable processes are essential to minimize the negative impacts of population growth. Bamboo is a viable alternative, because it is a product with potential properties for wood replacement. However, to achieve a good yield, it is necessary to identify suitable sites for the forestry of the species. This study aimed to forecast favorable areas for the best establishment and development of Bambusa vulgaris and Dendrocalamus giganteus in the Brazilian territory, functioning as a zoning for edaphoclimatic aptitude that can support decisionmaking. It was performed using georeferenced points of natural occurrence for the species, collected from a virtual herbarium. The environmental variables were extracted from the WorldClim, MOD16 and SoilGrids databases. In order to select the set of variables to be used for modeling each species, a correlation test was performed. From this result, ten environmental variables were used for each species, and non-correlated variables were chosen. A principal component analyses was applied to select the variables that better explained variability. For the distribution map, the Maxent algorithm was used, based on the principle of maximum entropy, evaluated by the AUC and TSS metrics. Bambusa vulgaris presents a large area for adaptation in the Brazilian territory, unlike the Dendrocalamus giganteus species, which showed a limited aptitude area.
Wood raw materials have showed scarcity in some Brazilian regions (Calegari et al. 2007, Biazus et al. 2010. Bamboo is a renewable resource of rapid growth and high yield that, in this context, could become a viable alternative to solve the limiting wood supply in these Brazilian regions (Calegari et al. 2007). However, to achieve a good yield, it is necessary to understand and identify suitable sites for the forestry of species of commercial interest (Nappo et al. 2005).
Bambusa vulgaris and Dendrocalamus giganteus are native from Asia. They have a wide variety of uses in construction, including houses and huts, small boats, fences, scaffolding, rustic furniture, handicrafts and cellulose (Schröder 2011). However, climate requirement data for these species are only general and superficial. In this context, climatic aptitude studies for the species introduction are essential.
Several factors define the plant growth; some are intrinsic to the species, others are related to the environment. The environmental factors are determinant for a good adaptability and development of a particular species in a site, with emphasis on luminosity, soil fertility, temperature and precipitation (Santos et al. 2011).
The predictive spatial distribution modeling of species is a technique that has been used to understand the climatic characteristics that determine its occurrence, generating information necessary for extrapolation to possible potential areas with its presence (Anderson et al. 2003). This technique relates species observations data to environmental predictor variables using statistical approaches, resulting in an ensemble of ecological requirements of the species (Guisan et al. 2017). It has been applied for multiple reasons, such as climate change distribution impacts (Wrege et al. 2017), contribution to priority conservation areas (Pompeu et al. 2018), biogeographic studies (Werneck et al. 2012) and as a tool for the extrapolation of possible potential areas for species cultivation, as in the study by Garcia et al. (2014).
Thus, observing the required use of renewable and sustainable materials to mitigate the negative impacts of population growth, in combination with the knowledge about climatic aptitude for the forestry of species of commercial interest, this study aimed to forecast favorable areas for the best establishment and development of Bambusa vulgaris and Dendrocalamus giganteus in the Brazilian territory, functioning as a zoning for edaphoclimatic aptitude that can support decision-making.

MATERIAL AND METHODS
The study was performed between 2016 and 2018, using two bamboo species: Bambusa vulgaris Schard. and Dendrocalamus giganteus Munro.
To run the model, it is necessary to get geographic coordinate points and environmental spatial database. For the modeling process, information are extracted for each point applied from the overlapped environmental variables, considering the limitation factors of each variable. The process is detailed below.
To compose the initial database, it was necessary to collect natural occurrence points of the species under study and edaphoclimatic layers. Initially, all georeferenced points of natural occurrence of the species used in this study were collected. The points were filtered, excluding inaccurate and doubtful data (Elith & Leathwick 2009), and, after the selection filter, the number of georeferenced points of natural occurrence were 17 for B. vulgaris and 12 for D. giganteus (Figure 1). These are satisfactory numbers, considering the study by Wisz et al. (2008), which addresses good results using Maxent for only 10 points. These georeferenced points were obtained from the herbarium database of the Global Biodiversity Information Facility (GBIF), a research vehicle widely used for studies with such purpose (Giannini et al. 2012).
Then, the collection of edaphoclimatic layers began. The first collections included 19 environmental layers related to precipitation and temperature of the Global Climate Data network, version 1.4 (Worldclim), representing interpolations of data observed from 1950 to 2000, with spatial resolution of 2.5 arc-min (~5 km) (Hijmans et al. 2005). Subsequently, the environmental layer of soil classification based on the FAO's World Reference Base was collected. This layer was obtained through the World Soil Information database (Soil Grids), a project that provides global forecast files for soil characteristics and universal classification, with a resolution of 1 km (Hengl et al. 2014). The soil layer is considered analyzing the similarity of classification between the origin sites of bamboo occurrence and the Brazilian soils through the universal classification, where the classification of soils is related with the soil properties. Environmental layers of potential evapotranspiration (PET) and actual evapotranspiration (AET) were also collected, with a spatial resolution of 1 km, obtained from the Moderate Resolution Imaging Spectro radiometer (MOD16) network, a project related to NASA studies. These data are widely used in some studies to estimate the evapotranspiration of the Earth's surface, using remote sensing data (Mu et al. 2011, Liu et al. 2013. After collecting all the mentioned layers, the water deficit (WD) was also inserted into the initial database. This layer was specifically created based on PET and AET, using the equation proposed by Pereira et al. (2002), where WD = PET -AET, using tools from the ArcGIS software, version 10.3. In this way, the initial database was formed, containing all the collected variables.
To avoid modeling conflicts, the spatial resolution of soil classification layers and those related to evapotranspiration and water deficit were changed to 2.5 arc-min (~5 km), thus being compatible with all variables in the initial database.
The initial choice of variables occurred as it follows: firstly, it was established that easily accessible environmental variables for different audiences were inserted in the set of variables that would be part of the modeling, being these: BIO1 (annual mean temperature) and BIO12 (annual precipitation, in mm); secondly, the water deficit was added as an environmental variable with a great impact on the survival factor for several agricultural and forest crops. In this way, based on the initial database, the variables BIO1, BIO12 and water deficit were established for the modeling of both species and treated, in this context, as primary variables.
Next, in order to guide the choice of other environmental variables that would be part of the modeling of each species under study, a statistical correlation analysis was performed among all climatic variables of the initial database. This statistical process was performed based on values extracted from each species occurrence points of all the involved edaphoclimatic variables. The test was performed using the Statistica software, version 12 (StatSoft 2013), applying the Pearson's correlation coefficient, generating a correlation table among all the environmental variables. For this study, it was established that values ≥ 0.8 (positive or negative) represented a high correlation. Finally, the collinear variables were removed, and the non-correlated variables were chosen to run the final model. In this way, it is guaranteed that the chosen variables are not redundant and that they represent more broadly factors that can be related to the occurrence of the species. This procedure was repeated separately for each species under study. The primary variables can be observed in Table 1, together with the variables selected through the correlation analysis for each referred species. At the end of this process, the data reliability is verified, in order to prevent biased information.
After the selection and definition of the set of edaphoclimatic variables used to model each species, a principal component analysis was performed using the Statistica software, version 12 (StatSoft 2013). The objective of the principal component analysis is to select the variables that most explain the variability inside the multidimensional dataset, in order to infer a Bio1: annual mean temperature; Bio2: mean diurnal range (monthly mean: maximum temperature -minimum temperature); Bio3: isothermality; Bio7: temperature annual range; Bio10: mean temperature of the warmest quarter; Bio12: annual precipitation; Bio16: precipitation of the wettest quarter; Bio17: precipitation of the driest quarter; Bio18: precipitation of the warmest quarter; Bio19: precipitation of the coldest quarter; WD: water deficit; SC: soil classification. factor of its importance, in relation to the occurrence and adaptation of the mapped species. The Maxent software, version 3.4.0 (Phillips et al. 2017) (or maximum entropy algorithm) was used to analyze the points of species occurrence in the suitability distribution. The software used the selected variables and georeferenced points of occurrence from the referred species. This method gives a greater weight to the points of species occurrence in suitability distribution and has benefits, in relation to other algorithms, because it requires only the presence of data from the species (Phillips et al. 2006). Other important benefit for the Maxent used in this study is the response curves outputs from this software, which show how the model prediction changes as each variable is varied, offering an important tool for assessing the biological plausibility of the model (Merow et al. 2013) and providing, for the present study, a supporting tool for cultivation decision-making using local available climate data. However, it is important to mention that it is better to apply the final map to choose cultivation areas where there is a higher percentage in adaptation suitability, because the final map integrates the set of edaphoclimatic variables.

Species
With the aid of this software, the points of species occurrence were crossed to the selected edaphoclimatic layers, extracting the existing information in each layer for each georeferenced point. With these extracted data, the software creates an information bank by the "machine learning" method, that is based on the species natural occurrence, in order to estimate the species adaptation suitability to the place of interest, based on the suitability of maximum entropy (Magalhães 2013). Using this information bank, Maxent creates two data sets called training and testing. Training data sets are used to train the software on the behavior of species occurrence related to edaphoclimatic variables, while the testing set is used to measure the effectiveness degree of the learned concept. Here, the cross-validation partition data was applied, randomly splitting the data into training and testing data (Franklin 2010), thus generating the final model that corresponds to the mean of 10 replications performed by the software.
The final model for each species was evaluated by two metrics: area under the curve (AUC; Barry & Elith 2006) and true skill statistic (TSS; Allouche et al. 2006), strongly applied in this type of species distribution study. The AUC metric ranges from 0 to 1, where 1 indicates a perfect discrimination, 0.5 a discrimination not better than a random assumption, and values lower than 0.5 imply a performance worse than the random assumption (Phillips et al. 2004, Giannini et al. 2012. TSS takes both omission and commission errors into account and ranges from -1 to +1, with +1 indicating a perfect agreement and 0 to -1 a random performance (Allouche et al. 2006). Here, were considered as successful final species models those presenting AUC > 0.8 and TSS > 0.6 (Gallien et al. 2012).
Through the map produced by Maxent and with the aid of the ArcMap software, version 10.3, a binary map was produced, in which the values of each pixel are represented by 0 or 1, indicating non-apt and apt areas, respectively, for the referred species. A logistic threshold of minimum training presence was used, thus covering all the areas that can be suitable for the cultivation of the studied species. Finally, with the aid of the ArcMap 10.3 software and using the raster calculation tool, a subdivision of areas considered apt in adaptation rates was extracted from the binary map.

RESULTS AND DISCUSSION
The potential area modeled for Bambusa vulgaris showed a good performance for AUC (0.92) and TSS (0.68). Through the generated map, the species suitability distribution throughout the Brazilian territory is observed ( Figure 2). Thus, policy makers can use the map to indicate areas for

cultivation of this species and farmers can check this map under their target lands.
Bambusa vulgaris presented an area with occurrence suitability (≥ 0.76) of approximately 600,000 km². Most of this area is in the northern region of the country, more precisely in the Amazon biome. However, part of this area also extends along the coastline of Brazil and stands out in the states of Bahia and São Paulo, where the Atlantic Forest biome is located. The 600,000 km² represent ten times the area of forest planted in the Brazilian territory (IBÁ 2017), it is an extensive area and, if there is interest in the expansion of the Bambusa vulgaris planted area in Brazil, it should begin with areas already anthropized, in regions with occurrence suitability ≥ 0.76.
According to the Köppen classification and Alvares et al. (2013), in the climatic zones of the Brazilian territory, it is observed that the species is more likely to occur in tropical climates without a dry season (Af). On the other hand, with a lower suitability, the species may also occur in tropical climates with a dry winter (Aw), monsoon climate (Am), oceanic climate without a dry season and with hot summer (Cfa), and oceanic climate without a dry season and warm summer (Cfb). This fact also corroborates Lobovikov (2007), which cites that the species is endowed with a great plasticity and capacity to develop in different areas, and still tolerates extreme conditions. The species is widely used throughout the country by small businesses, mainly as a source for energy generation. In addition, it is highlighted the large species plantation for pulp production destined to the manufacture of cement bags by the João Santos holding company, in Coelho Neto, Maranhão state, Brazil (Santi 2015), whose production site is mapped in our study with an adaptation suitability of up to 25 %, thus increasing the species expectations for areas that present greater probabilities. It is extremely important to know the contribution of each variable to the modeling, in order to visualize those that most contributed to the mapping of adaptation suitability. In this context, the contribution of each variable can be observed in Table 2.
The response curves of the four variables that most contributed to the modeling distribution of Bambusa vulgaris were annual precipitation (25.91 %), mean diurnal range (12.96 %), precipitation of the warmest quarter (11.94 %) and mean temperature of the warmest quarter (10.95 %) (Figure 3). Local climate data can be assessed with these response curves, helping in the decision-making process for the cultivation of this species. The annual precipitation curve ( Figure 3A) shows an exponential increase in the species adaptation suitability up to a precipitation of 3,000 mm, and then this suitability tends to stabilize. Thus, the importance of a high precipitation for the species is verified, explaining better the possible reason for the non-occurrence of the species in semiarid climate regions, where, according to Köppen, the annual precipitation might be lower than 800 mm.
On the other hand, the curve for mean diurnal range ( Figure 3B) may better explain the nonoccurrence of the species in places with a high range at some time of the year, such as in the South, Southeast and Center-West states (Alves & Barduco 2011, Souza et al. 2013, besides the highest occurrence suitability in the North region, where the diurnal range is low, in relation to other Brazilian regions. The response curve for the precipitation of the warmest quarter ( Figure 3C)  of almost 100 % for adapting Bambusa vulgaris in areas with values equal or greater than 2,000 mm. In this way, it is suggested a better species adaptation in places where the hot season is associated with a high precipitation, thus limiting the species adaptation in the northeastern regions. For the mean temperature of the warmest quarter ( Figure 3D), it is observed that values close to 38 ºC presented the best adaptation results for the species (80 %), and the decrease of this temperature reflects a significant fall in the adaptation suitability.
It is noteworthy that the variables water deficit and annual mean temperature did not contribute greatly to the distribution mapping of the species, in relation to other variables. Moreover, it is observed, by the contribution of the soil classification variable, that it is not a highly limiting factor for this species, what can further increase the adaptation suitability areas.
Dendrocalamus giganteus showed a limited aptitude area in the Brazilian territory (Figure 4). The statistics showed a good performance for AUC (0.94) and TSS (0.83). The adaptation area with the best adaptation suitability (≥ 0.76) accounts for 20,000 km².
As an exotic species, the presented result may reflect that some edaphoclimatic characteristics existing in its place of natural occurrence, not detected in the modeling, can be a limiting factor for the species establishment in places without the occurrence of such characteristics.
It can be observed that the mapped area is mainly comprised on the coast of the southeastern  region, more precisely in the northeast of the São Paulo state and in large part of the Rio de Janeiro state, where, according to the Köppen classification, predominate climates tropical with dry winter (Aw) and temperate humid with dry winter and warm summer (Cwb). These results corroborate those by Judziewicz et al. (1999), which mention that the species adapted well to the climate in the southeastern region, although this study does not demonstrate that the South region is suitable for the species occurrence, as mentioned by the authors.
The four variables that most contributed to the species distribution modeling were the precipitation of the coldest quarter (23.20 %), precipitation of the driest quarter (19.84 %), mean diurnal range (16.55 %) and water deficit (14.46 %) ( Table 3).
It is noteworthy that the variables annual mean temperature and annual precipitation were not strong contributors for mapping this species. However, the variable water deficit, which was initially highlighted as a possible limiting factor for several forest and agricultural crops, contributed significantly to the species distribution modeling.
The response curves of the four variables that most contributed to the modeling explain better the prediction result of the species occurrence, considereing that the following climates occur in these sites: tropical with dry winter (Aw) and humid temperate with dry winter and warm summer (Cwb), according to the zones described by Köppen ( Figure 5).
It was observed that the species tends to develop better in areas with a lower precipitation in colder ( Figure 5A) and dry seasons ( Figure 5B), considering that the increase of precipitation in these conditions is a limiting factor in the species adaptation suitability. Similar characteristics are showed by Tectona grandis, which has a betterquality development when subjected to periods of mild temperatures associated with dry periods, from three to five months (Remade 2004).
Analyzing the response curve for the mean diurnal range ( Figure 5C), it is notable that the species shows a certain adaptation sensitivity to its increase. The best species adaptation (≈ 65 %) is observed in the temperature of 2-5 ºC, and a diurnal value above this range may decrease the species adaptation suitability.
In this way, these results may explain the lack of species adaptation in other Brazilian regions, which can show climate ranges with high temperatures throughout the year or summer (North and Northeast); high mean diurnal ranges (Center-West) (Alves & Barduco 2011, Souza et al. 2013; and places with climates without dry seasons (South).
For the water deficit, it was observed a high range of adaptation suitability related to water deficit (0-350 mm), where values above 150 mm gradually decrease the adaptation suitability, and values close to 100 mm show a better adaptation for the species. This species can support higher values of water deficit (around 350 mm), showing that it is apt to occur in regions where there is a dry season.
Although our model shows a high accuracy for both species, it is important to highlight that the environmental suitability was predicted based on the natural occurrence database of these species, showing the optimum sites for its occurrence through environmental database. However, these species could not occur beyond their natural location, due to ecological issues (e.g.: geographical barriers, competition and dispersion), rather than environmental constraints. This may cause a higher level of suitability constraint on the adaptation map, since the modeling is based on environmental data similarity. Other point is that our database shows few geo-referenced points. So, these species could occur in locations where there is no global database available. Thus, the locations with high suitability show a high environmental similarity with the natural occurrence of bamboo, being highly apt for cultivation. Other locations are not necessarily unfit, but should be analyzed with caution, and growth monitoring should be carried out in these locations, before making commercial investments.
There are scarce or superficial studies of these species related to the adaptation of climatic conditions, what hinders the discussion of the presented results. However, this study supports the decision-making for the best establishment and development of B. vulgaris and D. giganteus in Brazil.
CONCLUSIONS 1. Bambusa vulgaris shows greater aptitude areas for cultivation in the Brazilian territory, in relation to Dendrocalamus giganteus, and B. vulgaris is more suitable to occur in northern, northeastern and coastal regions of the country; 2. For the Center-West region and the Brazilian Savannah areas, the studied species show little adaptation, being necessary monitoring their growth before considering cultivation investments; 3. B. vulgaris is more apt for cultivation in regions where there are higher values of precipitation and temperature, while D. giganteus develops better in regions where there are dry seasons.