The Effects of Climate Change on the Development of Tree Plantations for Biodiesel Production in China

Biodiesel produced from woody oil plants is a promising form of renewable energy but a combination of tree plantations’ long cultivation time and rapid climate change may put large-scale production at risk. If plantations are located in future-unsuitable places, plantations may fail or yield may be poor, then significant financial, labor, and land resources invested in planting programs will be wasted. Incorporating climate change information into the planning and management of forest-based biodiesel production therefore can increase its chances of success. However, species distribution models, the main tool used to predict the influence of future climate–species distribution modeling, often contain considerable uncertainties. In this study we evaluated how these uncertainties could affect the assessment of climate suitability of the long-term development plans for forest-based biodiesel in China by using Sapindus mukorossi Gaertn as an example. The results showed that only between 59% and 75% of the planned growing areas were projected suitable habitats for the species, depending on the set-up of simulation. Our results showed the necessity for explicitly addressing the uncertainty of species distribution modeling when using it to inform forest-based bioenergy planning. We also recommend the growing area specified in China’s national development plan be modified to lower the risk associated with climate change.


Introduction
Rapid deployment of renewable energy has been proposed as a way to meet the obligation of greenhouse gas reduction stated in the Paris Agreement [1,2].Forest-based bioenergy is a promising form of renewable energy because it is capable of providing multiple co-benefits [3].Well planned and managed forest-based bioenergy production can help to enhance energy and food security, create job opportunities, and reduce environmental degradation [4][5][6][7].Forest-based bioenergy has, therefore, been given increased attention in many countries and regions, such as the United States [8], Europe [9], and China [10].
Among the different forms of forest-based bioenergy, biodiesel produced from woody oil plants is especially attractive because it can serve as a substitute for fossil oil [11].Biodiesel can be used directly in modified diesel engines or blended with conventional diesel [12].Development of biodiesel production based on woody oil plants at a large scale needs careful planning.Growers are generally Forests 2017, 8, 207 2 of 16 advised to match species and cultivars to local climates and site conditions but the possible impact of mid-to long-term climate change was seldom considered in the development plan [10].Nevertheless, significant climate changes including warming and shifted precipitation patterns have been projected for the 21st century [13].Woody oil plants normally need 6 to 20 years to start producing fruit and will stay productive for several decades, e.g., the species Xanthoceras sorbifolia Bunge reaches full-production age around 20 years old and has a productive life of over 100 years [14].Jatropha curcas L. can be productive for more than 30 years [15].This long production time means any new energy tree plantations planted according to a long-term development plan will be subjected to the impact of rapid climate change.If plantations are located in places where future climate is unsuitable for the planted species, plantations may fail or the yield may be poor.As a consequence, the significant financial, labor, and land resources invested in the planting programs will either be wasted or not have the expected biodiesel yield.
A few studies have already shown the potential impact of climate change on woody oil plants.For example, the area suitable for growing Millettia pinnata (L.) Panigrahi would extend southward into the subtropics in Australia [16].Global warming was predicted to affect suitable areas for three woody oil plants (Sapium sebiferum (L.) Roxb., Vernicia fordii (Hemsl.),X. sorbifolia) in China negatively by 2080 under the IPCC scenario A2a [17].Seed yield of Jatropha was predicted to decrease in the Sahel, Eastern Brazil, and Northern Australia due to climate changes [18].While these studies all provided useful information, the incorporation of climate change information into the policy making process of forest-based bioenergy needs more attention on an important issue, i.e., the uncertainty in projections of suitable habitats for woody oil plants.Species distribution modeling (SDM) is the primary tool used to predict how species' distribution ranges will change in the wake of future global climate change [19].It is well known that variability of data quality, modeling methods, general circulation models (GCMs), emission scenarios, thresholds, and dispersal ability can influence the magnitude, direction, and rates of change in projections of suitable climate space of the modeled organism and influence the outcome of targeted management [20][21][22][23].These sources of uncertainty and their influences on projections have not been addressed adequately.A relevant issue is how to incorporate the variability of projections into the planning of biodiesel production form woody oil plants.Currently, this issue was handled in various ways in forest planning, from rejecting the usefulness of the SDMs for informing policy planning [24], to finding consensus projections [25], and to adopting adaptive forest management programs that consider the range of possible trends and uncertainties [26,27].Studies on these issues are urgently needed in order to prepare forest-based bioenergy production for future climate change.
China holds great potential to develop forest-based biodiesel.A total of 154 woody oil plant species with fruits or seeds containing more than 40% fat content, have been identified in China thus far [28].Along with the abundant plant species, China has more than 43.75 million ha of marginal lands that can potentially be used for growing energy tree plantations [29].The State Forestry Administration (SFA) of China issued the Mid-term Development Plan for Forest-based Bioenergy to guide the development of forest-based bioenergy in China between 2012 and 2020 [30].In this plan, the total area allocated for woody oil plant plantations would be increased from 1.35 million ha to 4.22 million ha and the potential yield of biodiesel would reach 580 million tons of standard coal equivalents by 2020.In order to meet these targets, the SFA planned to construct large production bases containing tree plantations and processing facilities in specified regions.An examination of the climate suitability of China's development plan will not only inform policy makers in China but also can provide lessons for other countries where large-scale production of forest-based biodiesel is under consideration.
Six woody oil plant species including X. sorbifolia, Elaeis guineensis Jacq., Pistacia chinensis Bunge, J. curcas, Swida wilsoniana (Wanger.)Sojak, and Sapindus mukorossi Gaertn have been chosen for large-scale planting in China's development plan.In this study we focused on S. mukorossi (soapberry), for which we have conducted extensive field surveys on its distribution in China.We wanted to find out how predictions of the distribution range of S. mukorossi varied when different aspects of uncertainty (GCMs, emission scenarios, and dispersal ability) of SDM were addressed.We also wanted to demonstrate how the climate suitability of China's development plan would be affected by the variability of projections.

Compilation of Occurrence Data
Occurrences of S. mukorossi were compiled from various floras, books, and herbarium records.A list of sources used in this study can be found in the supplementary materials (Table S1).The compiled occurrence data can be divided into two types: (1) records with names of places and geodetic coordinates where the species were identified; and (2) records with names of places but no geodetic coordinates.For the second type of records, the center coordinates of the smallest administrative units corresponding to the names of places were used as the coordinates of the record.Occurrence records of species often contain spatial bias that can affect the accuracy of model predictions [31,32].In this study a spatial filtering approach was used to reduce spatial bias by merging all occurrence records in the same grid (size = 10 km) as one occurrence.In total, 734 out of 924 occurrence records were kept after filtering.
We developed an independent dataset of occurrences using citizen science data and data we collected in field surveys.The citizen science data were extracted from Plant Photo Bank of China (PPBC), a web-based system to collect, collate and preserve plant photographs, slides, and digital images taken by researchers as well as botany enthusiasts.All images included in the system have been validated by botanists from the Institute of Botany, the Chinese Academy of Sciences [33].Data in PPBC included observations made from 2008 to the present.We conducted field surveys in the main distribution region of S. mukorossi in China from July 2011 to November 2013.Coordinates of places where trees of S. mukorossi were found in natural habitats were recorded using handheld GPS units.A total of 76 and 50 records were obtained from PPBC and our field surveys, respectively.

Climatic and Environmental Inputs
Climatic variables including both the base and projected climate data sets were acquired from the online database, WorldClim [34].The base set is the representative climate from 1950 to 2000.The projected set are projections for future climate between 2041 and 2060 (2050s) from 18 general circulation models (GCMs) included in Coupled model intercomparison project phase 5 (CMIP5).The names of the GCMs used in this study and their producers can be found in the supplementary files (Table S2).The 2050s was selected as the targeted period because S. mukorossi tree plantations planted according to China's development plan will reach the stage of full production by that time.Two representative concentration pathways, RCP2.6 and RCP8.5, were considered in this study.RCP8.5 has the highest radiative forcing values (8.5 W/m 2 ) while RCP2.6 has the lowest radiative forcing values (2.6 W/m 2 ) by 2100.Together they can provide a range of future climate change scenarios.Both base and projected data sets contained 19 bioclimatic variables.
Four biophysical variables including altitudes, slopes, aspects, and soil classifications were extracted from the Harmonized World Soil Database [35].Existing studies showed they had significant effects on predicting distributions of plant species at the regional and continental scale [36,37].The selected bioclimatic and biophysical variables were converted to Raster ASCII grid format with a resolution of 5-arc min.

Model Simulation
The MaxEnt version 3.3.3ksoftware [38] was used to predict the distribution ranges of S. mukorossi.MaxEnt requires only species occurrence data and environmental variables and has performed well in recent comparative studies [39,40].The extent of the study area encompassed 16 provinces and one municipality in eastern China.This extent was used for modeling the potential distribution of S. mukorossi because all occurrence records fell within this region.The 10 provinces and one municipality where large-scale plantations have been planned by the SFA were also located in this study area (Figure 1).The SFA planned to add 0.25 million ha of new plantations of S. mukorossi by 2020 [30].have been planned by the SFA were also located in this study area (Figure 1).The SFA planned to add 0.25 million ha of new plantations of S. mukorossi by 2020 [30].To predict the current potential distribution of S. mukorossi, the compiled occurrence data were randomly split into 70% for training model and 30% for validating model results.To reduce the collinearity of variables, pre-model runs were conducted using all variables, and highly crosscorrelated variables were removed following the method used by Garcia et al. [41].
The area under curve (AUC) was used to assess the accuracy of the model.AUC values > 0.7 are considered to produce "useful" model outputs with good discriminatory power [42,43].The continuous predictive values were converted into binary data sets by using threshold selection techniques.Two methods were tested, including maximizing the sum of sensitivity and specificity (MSS) and equalizing the sensitivity and specificity (ESS) [44].Predicted probabilities greater than the threshold value were assigned the value of one, representing moderate to high suitability of habitats (hereafter suitable habitats).Prediction values less than the threshold value were assigned the value of zero, representing low to unsuitable habitats (hereafter unsuitable habitats).The percentages of independent observations that have been correctly predicted were used to judge the accuracy of potential current distribution maps and select the thresholding method used for further analysis.
The future distribution of S. mukorossi was predicted using all occurrence records (compiled records + independent data) and future climate projections.Because relying on only one or few arbitrarily-selected climate projections can increase the likelihood of producing biased projections [45], simulations were conducted in two different ways: (1) averaged bioclimatic variables retrieved from the outputs of 18 GCMs were used in simulations; and (2) the output of each GCM was used in separate simulations and the results of multiple individual projections were combined to generate a consensus projection map [25].Two assumptions on species dispersal, i.e., full dispersal and no dispersal, were used [46].Eight projections were produced at the end (two ways of using outputs of To predict the current potential distribution of S. mukorossi, the compiled occurrence data were randomly split into 70% for training model and 30% for validating model results.To reduce the collinearity of variables, pre-model runs were conducted using all variables, and highly cross-correlated variables were removed following the method used by Garcia et al. [41]. The area under curve (AUC) was used to assess the accuracy of the model.AUC values > 0.7 are considered to produce "useful" model outputs with good discriminatory power [42,43].The continuous predictive values were converted into binary data sets by using threshold selection techniques.Two methods were tested, including maximizing the sum of sensitivity and specificity (MSS) and equalizing the sensitivity and specificity (ESS) [44].Predicted probabilities greater than the threshold value were assigned the value of one, representing moderate to high suitability of habitats (hereafter suitable habitats).Prediction values less than the threshold value were assigned the value of zero, representing low to unsuitable habitats (hereafter unsuitable habitats).The percentages of independent observations that have been correctly predicted were used to judge the accuracy of potential current distribution maps and select the thresholding method used for further analysis.
The future distribution of S. mukorossi was predicted using all occurrence records (compiled records + independent data) and future climate projections.Because relying on only one or few arbitrarily-selected climate projections can increase the likelihood of producing biased projections [45], simulations were conducted in two different ways: (1) averaged bioclimatic variables retrieved from the outputs of 18 GCMs were used in simulations; and (2) the output of each GCM was used in separate simulations and the results of multiple individual projections were combined to generate a consensus projection map [25].Two assumptions on species dispersal, i.e., full dispersal and no dispersal, were used [46].Eight projections were produced at the end (two ways of using outputs of GCMs × two climate scenarios × two assumptions of dispersal).For projections using averaged GCM outputs, the projected probability maps were converted to binary distribution maps by using the selected thresholding method.For projections using outputs of individual GCMs, the consensus distribution maps were generated using a majority voting approach.After converting each projected probability map into the binary suitable/unsuitable map, the grid cell was assigned to the suitable (or unsuitable) category if more than 50% of the individual projections were the same.
All simulations were implemented using a cross-validation approach included in MaxEnt with ten replicates.The entire modeling process is depicted in the following flowchart (Figure 2).GCMs × two climate scenarios × two assumptions of dispersal).For projections using averaged GCM outputs, the projected probability maps were converted to binary distribution maps by using the selected thresholding method.For projections using outputs of individual GCMs, the consensus distribution maps were generated using a majority voting approach.After converting each projected probability map into the binary suitable/unsuitable map, the grid cell was assigned to the suitable (or unsuitable) category if more than 50% of the individual projections were the same.All simulations were implemented using a cross-validation approach included in MaxEnt with ten replicates.The entire modeling process is depicted in the following flowchart (Figure 2).

Analysis of Model Outputs
In order to assess the impact of climate change on the distribution of S. mukorossi, the predictions of current potential distribution were compared with the predicted future distribution maps (2050s).The loss and gain of areas suitable for the species were determined through comparisons.
To determine the difference among future projections caused by different set-ups of the model, the predicted distribution maps were compared using Kappa coefficient of agreement, an indicator used for comparing two maps [47].The values of Kappa range from 1 to −1 with −1 representing no overlap between two maps and 1 representing perfect overlap between the two maps.The value zero represents the special case where the agreement is equal to the agreement that can be expected by chance [47].
Finally, the predicted suitable areas from various models were compared with the planned growing areas of S. mukorossi specified in the SFA's development plan.Percentages of overlapped areas to the total planned growing areas were calculated to indicate the agreement between projected habitats and planned growing areas.

Selection of Model Variables and the Modeling Accuracy
Based on their contributions to model output calculated in the pre-model run and the correlation values, six bioclimatic variables were used in MaxEnt, i.e., minimum temperature of the coldest month, isothermality, temperature seasonality, maximum temperature in the warmest month, annual precipitation, and precipitation seasonality.All four biophysical variables (previously mentioned) were kept for further analysis.
The values of AUC of all model projections were higher than 0.7 (Table 1), indicating acceptable accuracy of modeling.

Analysis of Model Outputs
In order to assess the impact of climate change on the distribution of S. mukorossi, the predictions of current potential distribution were compared with the predicted future distribution maps (2050s).The loss and gain of areas suitable for the species were determined through comparisons.
To determine the difference among future projections caused by different set-ups of the model, the predicted distribution maps were compared using Kappa coefficient of agreement, an indicator used for comparing two maps [47].The values of Kappa range from 1 to −1 with −1 representing no overlap between two maps and 1 representing perfect overlap between the two maps.The value zero represents the special case where the agreement is equal to the agreement that can be expected by chance [47].
Finally, the predicted suitable areas from various models were compared with the planned growing areas of S. mukorossi specified in the SFA's development plan.Percentages of overlapped areas to the total planned growing areas were calculated to indicate the agreement between projected habitats and planned growing areas.

Selection of Model Variables and the Modeling Accuracy
Based on their contributions to model output calculated in the pre-model run and the correlation values, six bioclimatic variables were used in MaxEnt, i.e., minimum temperature of the coldest month, isothermality, temperature seasonality, maximum temperature in the warmest month, annual precipitation, and precipitation seasonality.All four biophysical variables (previously mentioned) were kept for further analysis.
The values of AUC of all model projections were higher than 0.7 (Table 1), indicating acceptable accuracy of modeling.

Current Potential Distributions of S. mukorossi
Using the independent data as the validation data, we obtained an accuracy rate of 73% for the binary distribution map produced using MSS method and 61% for that produced using ESS.We kept the map produced using MSS as the current potential distribution map of S. mukorossi and used MSS as the thresholding method in future projections.
The current potential distribution of S. mukorossi displayed a southeast-to-northwest gradient (Figure 3a).The potential habitats for the species were most abundant in southeast China, including Guangdong Province, Guangxi Province, Fujian Province, and Zhejiang Province.The distribution of the species decreased when moving north and west.Sporadic distribution occurred in Sichuan province, Shaanxi Province, Henan Province, Shandong Province and Chongqing.

Current Potential Distributions of S. mukorossi
Using the independent data as the validation data, we obtained an accuracy rate of 73% for the binary distribution map produced using MSS method and 61% for that produced using ESS.We kept the map produced using MSS as the current potential distribution map of S. mukorossi and used MSS as the thresholding method in future projections.
The current potential distribution of S. mukorossi displayed a southeast-to-northwest gradient (Figure 3a).The potential habitats for the species were most abundant in southeast China, including Guangdong Province, Guangxi Province, Fujian Province, and Zhejiang Province.The distribution of the species decreased when moving north and west.Sporadic distribution occurred in Sichuan province, Shaanxi Province, Henan Province, Shandong Province and Chongqing.

Future Distribution of S. mukorossi
The projections showed that, under the assumption of full dispersal, the distribution range of S. mukorossi would expand under both climate scenarios (Figure 3b-e).The largest expansion, 33.6%, was projected for the RCP8.5 scenario by using the consensus projection method (Figure 3e).The smallest expansion, 14.5%, was projected for the RCP 2.6 scenario by using the average GCMs method (Figure 3b).The expansion was the net result of loss and gain of suitable habitats (Table 2).Under the assumption of no dispersal, the distribution range of S. mukorossi would shrink under both climate scenarios (Figure 3f-i).The largest shrinkage, 3.8%, was predicted under the RCP 2.6 scenario using the average GCMs method (Figure 3f).
In order to reveal the difference in projection results, we further compared the eight projected future distribution maps.Agreement among these maps varied between 68% and 99% (Table 3).The projection with the assumption of full dispersal, high emission scenarios, and average GCMs method in general had low agreement with other projections.

Comparison between the SFA's Plan and Future Distributions
In the SFA's plan, 10 provinces and one municipality were listed as areas for developing plantations of S. mukorossi.We found that the overlap rates among the planned growing areas and projected suitable habitats for the species varied between 59% and 75% (Table 4).Provinces in southeast China had high overlap rates under all modeling set-ups.The province at the northwest corner of the planned growing area, i.e., Shaanxi Province, had the lowest overlap rate.In general, the overlap rates among the planned growing areas and projected suitable areas were higher under the assumption of full dispersal than those under the assumption of no dispersal if other conditions were the same.Similarly the overlap rates were higher under the high emission scenario under the full dispersal scenario.Inside each province, the overlapped area varied spatially (Figure 4a-h).southeast China had high overlap rates under all modeling set-ups.The province at the northwest corner of the planned growing area, i.e., Shaanxi Province, had the lowest overlap rate.In general, the overlap rates among the planned growing areas and projected suitable areas were higher under the assumption of full dispersal than those under the assumption of no dispersal if other conditions were the same.Similarly the overlap rates were higher under the high emission scenario under the full dispersal scenario.Inside each province, the overlapped area varied spatially (Figure 4a-h).

Changes of the Distribution Range of S. mukorossi
Our result showed that, under the full dispersal scenario, habitats of S. mukorossi would expand northwest.The northwest shift found in this study was in agreement with the general observation that distributions of plant species would shift towards the poles under global warming [48,49].Our results were similar to findings of several studies on the effects of climate change on distribution of woody oil plants, such as Sapium sebiferum in the United States [50], Olea europaea L. in the Mediterranean basin [51,52], and Hevea brasiliensis (Willd.ex A. Juss.)Muell.Arg. in China [53] and in India [54].
If no dispersal was assumed, the net change of distribution range would be negative.This may be due to the biology of S. mukorossi, which likes full sun and can tolerate drought but not waterlogging [55,56].This was affirmed by analyzing the contribution of bioclimatic variables to the

Changes of the Distribution Range of S. mukorossi
Our result showed that, under the full dispersal scenario, habitats of S. mukorossi would expand northwest.The northwest shift found in this study was in agreement with the general observation that distributions of plant species would shift towards the poles under global warming [48,49].Our results were similar to findings of several studies on the effects of climate change on distribution of woody oil plants, such as Sapium sebiferum in the United States [50], Olea europaea L. in the Mediterranean basin [51,52], and Hevea brasiliensis (Willd.ex A. Juss.)Muell.Arg. in China [53] and in India [54].
If no dispersal was assumed, the net change of distribution range would be negative.This may be due to the biology of S. mukorossi, which likes full sun and can tolerate drought but not waterlogging [55,56].This was affirmed by analyzing the contribution of bioclimatic variables to the distribution.Among all bioclimatic variables, minimum temperature of coldest month had the largest (18-27%) and annual precipitation the second largest (12-22%) contribution to the distribution of S. mukorossi.The climate of our study area was projected to become warmer and more humid by 2050s.The rise of temperature would make more regions suitable for the species and the distribution range of the species would expand.Our projections showed that the expansion would mainly occur in east Sichuan, west Chongqing, southeast Henan, and North Anhui, mostly mountainous areas.In the wild, S. mukorossi were often found to form large tree groves on slopes [57].The temperature rise in those mountainous areas would, therefore, add new habitats to the species.On the other hand, increased precipitation may make some current habitats no longer suitable, especially in places where soil drainage is an issue.Poor growth, diseases and even death has been reported for trees of S. mukorossi planted on sites with poor drainage or waterlogging [58].

Uncertainties in Projections
Paired comparisons of eight future distribution maps had results ranging from nearly identical (99%) to largely the same (68%) among these projections, demonstrating the uncertainty associated with using the SDM to project climate change impact on tree species.The divergent projections confirmed the need to pay more attention to uncertainty associated with species distribution models [26,59].Underestimation or overestimation of climate change impact could occur if the uncertainties are not explicitly addressed.For example, if we only modeled the future distribution of S. mukorossi under the full-dispersal assumption, high emission scenario, using the consensus method, half of Chongqing would be considered suitable for growing plantations of the species.However, when the other sources of uncertainty are considered, we would exclude the region from the growing plan.
The divergence also showed the need to accommodate uncertainty when using SDM in forest management.We have paid particular attention to the source of uncertainty in this study.Unlike other studies, we did not choose the entire country as our study area [17,60] because MaxEnt results were significantly affected by the landscape used for the background sample [40,61].We used spatial filtering to address the spatial bias associated with the occurrence records compiled from literature and herbarium specimens [62].We used an independent data set to help tune the model and select the thresholding method to transform continuous distribution of species into binary predictions of presence/absence, which is often the least explored source of uncertainty in SDMs [63].Also, the use of data collected from well-planned surveys will reduce the uncertainty in projections significantly [64].To address the uncertainties associated with the climate model predictions [65] and emission scenarios [66], we tried two different ways to use GCM output in MaxEnt and two emission scenarios.Our results clearly show the infeasibility of reaching one single projection with high certainty at this stage.The range of possible trends must be accommodated in the planning and policy-making process.

Implications for Forest-Based Bioenergy Planning
Predictions of potential future species distributions normally contain considerable uncertainties.So the question is how forestry administrations can be informed by these divergent results.First of all, as pointed out by Wiens et al. [59], ignoring the future is not an option for management practices and decisions despite the uncertainty associated with SDM.Secondly, the uncertainty could be incorporated into the decision-making process by following the principle of "prepare for the worst" and accommodating the range of possible trends and uncertainties [26,27].Using China as an example, the current SFA's development plan for forest-based bioenergy ignored future climate change entirely and thus carries a huge risk.The results showed that in some regions, i.e., Shaanxi Province, Guizhou Province, and Chongqing, percentages of predicted suitable areas were consistently low using different modeling set-ups.Four provinces including Guangdong, Guangxi, Fujian, and Jiangxi were consistently predicted as suitable for growing S. mukorossi.Other provinces fell in between these two groups.Under the principle of "prepare for the worst", it would be prudent to reconsider the plan to develop plantations of S. mukorossi in Shaanxi, Guizhou, and Chongqing as these regions were consistently predicted to be less suitable to the species than other provinces under future climate.More consideration should be given to the four provinces that were predicted to be suitable for the species.In the remaining provinces, plantations of S. mukorossi should be developed in predicted suitable regions and with climate adaptation measures in place, such as equipping sites with well-designed drainage systems.

Limitations of the Current Study
Although we tried to address the multi-source uncertainty associated with SDMs and explored a way to incorporate uncertainty into forest-based biodiesel planning and management in China, we did not exhaust all sources of uncertainties.For example, some studies have reported the species distribution model itself was a factor contributing to uncertainty [21,22,63].We did not compare the MaxEnt model to other species distribution models.We also did not quantify the relative contribution of different sources of uncertainty, which has been proposed as a strategy to improve robustness of planning and decision-making [21,22].We also did not consider factors that can affect the realization of projected ranges of plans on the ground [67], e.g., the competition for lands with other uses and markets for biodiesel.However, since our main purpose was to provide a comparison between planned growing areas specified in the government plan and projected suitable habitats for S. mukorossi, and the comparison results clearly showed the divergence, the method used in this study could satisfy that purpose.Our results clearly demonstrated that China's current development plan carries risks under the future climate and needs amendments.

Conclusions
Taking climate into consideration in forest planning and management has long been proposed by researchers as well as forestry organizations [68][69][70][71].However, the proposal was rarely applied in the field, especially in less-developed countries.Through this study, we developed the procedure for using species distribution modeling (SDM) to assess the national development plan for producing forest-based biodiesel in China.We paid special attention to sources of model uncertainties, including the data quality, GCMs, emission scenarios, thresholding method and dispersal ability.Using S. mukorossi as an example, we showed that the mismatch between the future likely distribution of the species and the planned growing areas can be as high as 41% of the total area of the planned growing areas.Although projections were divergent, the results still provided useful information.They showed that the current development plan, which does not consider climate change, carried possible risks such as failure of plantations and poor return on investment, under the future climate.This result can serve as a wake-up call for countries where the development of large-scale forest-based bioenergy is planned.Our study also showed that in order to use SDM to incorporate climate change information into forest management its uncertainties must be explicitly addressed to avoid misinformation.

Supplementary Materials:
The following are available online at http://www.mdpi.com/1999-4907/8/6/207/s1,Table S1: Sources of species occurrence records; Table S2: General Circulation Models (GCMs) used in this study; Dataset S1: Sapindus mukorossi presence data and independent testing data.The "Sapindus presence data.csv"file includes the location of S. mukorossi from the online database and published books.The "TestData.csv"file includes 127 locations of S. mukorossi from field surveys (the former 47 records in the table) and the Database of Plant Photo Bank of China (PPBC) (the remaining 80 records), which are independent to the records of "Sapindus presence data".

Figure 1 .
Figure 1.Map of the study area.The locations of Sapindus mukorossi observances are shown as dots.

Figure 1 .
Figure 1.Map of the study area.The locations of Sapindus mukorossi observances are shown as dots.

Figure 2 .
Figure 2. The modeling process for generating future potential distribution of Sapindus mukorossi in China (based on MaxEnt).

Figure 2 .
Figure 2. The modeling process for generating future potential distribution of Sapindus mukorossi in China (based on MaxEnt).

Table 1 .
The AUC scores for MaxEnt output, in various modeled projections of Sapindus mukorossi distribution in China.

Table 1 .
The AUC scores for MaxEnt output, in various modeled projections of Sapindus mukorossi distribution in China.

Table 2 .
Predicted range shifts of suitable habitats in 2050s (2041-2060) relative to the current potential distributionfor Sapindus mukorossi in China.

Table 3 .
Values of Kappa coefficients between different projection maps of future potential distribution of Sapindus mukorossi in China.

Table 4 .
Overlaps between the projected suitable habitats and the areas of provinces selected for growing Sapindus mukorossi in China.