Introduction

Climate change has been documented at broad resolution scales since the mid XX century. Since then, huge quantities of resources have been invested in the understanding of the phenomenon and its possible consequences on species distribution (^{Martínez-Meyer, 2005}). Paleoecological evidence indicates that the distribution area of some species has contracted or expanded on various occasions since the last ice age; therefore, it is assumed that in the future there could be a similar phenomenon. The decrease or increase in the geographic distribution ranges of the species are in function of their adaptability to the climate change. Such adaptability is related to genetic variation; the species with less variation will be more susceptible to extinction (^{Hamrick, 2004}).

Some studies recommend that species sensitive to climate change must be examined closely in order to define preservation strategies ex situ (^{Sáenz-Romero, Rehfeldt, Duval, & Lindig-Cisneros, 2012}; ^{Viveros- Viveros et al., 2009}). An alternative study to evaluate the impact of climate change on the geographic space of the species is through the use of species distribution models or ecological niche models (^{Jeschke & Strayer, 2008}; ^{Peterson, 2006}).

The species distribution modeling is done with two techniques: mechanistic and correlative. The former is based on physiological experiments of a particular species, whereas the latter links the environmental variables with the areas of species occurrence (^{Anderson & Raza, 2010}). The generated niche model is applied to a geographic space in order to identify ideal areas for species distribution. This type of modeling utilizes predicted future climates based on a combination of general atmosphere-ocean circulation models and reports of various scenarios. Such models pose different hypotheses on the effect of oceans and the atmosphere on climate; the effect is characterized by different supposed futures of socioeconomic and technological development (^{Real, Luz Márquez, Olivero, & Estrada, 2010}).

An important aspect in modeling the species distribution is delimiting the area of study prior to modeling, especially when invasive species or climate change are investigated, so that the predictive models have a high transferability (^{Anderson & Raza, 2010}). In the case of this study, the state of Michoacán was selected due to it being one of the entities with greatest plant diversity in Mexico (^{Cué-Bär, Villaseñor, Arredondo-Amezcua, Cornejo-Tenorio, & Ibarra- Manríquez, 2006}; ^{Villaseñor, 2003}). Among them, the Pinaceae family represents 27 % of the total area of study and at least 15 species belong to the Pinus genus (^{Cué et al., 2006}; ^{Gutiérrez & Duivenvoorden, 2010}). The importance of these species lies in their instrumental and intrinsic value. The former refers to goods such as food, fuel, nutrient recycling, water filtration, and contamination adsorption, and it can acquire greater economic value if the possibilities of biotechnology in genetic engineering are considered in the future. With regard to the intrinsic value, the decrease in biological diversity represents the loss of genetic viability of the species and ecosystems, for example, the loss of tree species and, in turn, the loss of habitat for other species such as Danaus plexippus L. (monarch butterfly) (^{Meléndez, 2010}; ^{Sáenz-Romero et al., 2012}; ^{Villaseñor, 2005}). In Mexico, there is no adequate system to designate protection areas for Pinus species (56 species) sensitive to climate change, given that less than 10 % of the distribution area is designated as a natural protected area (^{Gutiérrez & Duivenvoorden, 2010}).

This paper had the objective of evaluating the impact of climate change on the potential distribution of 12 Pinaceae species (10 pines and two fir trees) in Michoacán, with the purpose of identifying the species that will be most affected.

Materials and methods

Area of study

The state of Michoacán (Figure 1) has an approximate surface area of 58,864 km². Altitudes vary from sea level up to 3,840 m (Tancitaro volcano); the average altitude is 1,343 m. The state is located inside three physiographic regions: the Balsas Depression, the Pacific Coast and the Trans-Mexican Volcanic Belt. The main types of rocks are basic and acid extrusive igneous (^{Instituto Nacional de Estadística y Geografía e Información [INEGI], 1988}). The most prominent types of soil are Vertisol, Luvisol, Latosol, Andosol, Regosol and Phaeozem (^{INEGI, 1988}). The electric conductivity of the soils is less than 0.22 dS·m^-1 and the content of organic matter varies from 1.4 to 3.4 %; the greatest soil ratio has a slightly to moderate acidic pH (^{Cruz-Cárdenas et al., 2014}). The climate types are warm sub-humid (Aw₁, Aw₀) and temperate sub-humid (C(w₁), C(w₂)) (^{García, 1990}); the minimum and maximum temperature averages are 13.8 and 29.4 ºC, respectively, with an average precipitation of 850 mm. The main plant communities found in Michoacán are: coniferous forests, oak forests, mountainous cloud forests, deciduous tropical forests, subtropical thickets, thorn forests and tropical semi-deciduous forests (^{Villaseñor, 2005}).

Figure 1 Area of study of the potential distribution of 12 species of Pinaceae in the state of Michoacán. 

Results and discussion

The climate variables selected with multinomial regression were: maximum temperature of the warmest month (bio5), minimum temperature of the coldest month (bio6), annual temperature range (bio7), and precipitation of the driest quarter (bio17). Based on the X² test, the Akaike information criterion for the model with 19 variables was not significantly different from the model with selected variables (P = 0.999). In regard to the model with the seven soil variables and the model with the selected soil variables (pH and K), there was also no significant difference (P = 0.395) between them. Therefore, in the species distribution modeling, four climate variables, two soil variables, and the digital elevation model were utilized. The importance of variable selection lies in reducing the co-linearity and thus preventing the risk of over adjusting the models. The selected climate variables in the temperature and precipitation groups are in accordance with what was proposed by ^{Barbet-Massin and Jetz (2014)}, with which models consistent in space and time are generated, mainly with the variables in the temperature group.

Table 2 shows the values of the area under the curve of the species distribution models. In general, the values of the area in all the scenarios are greater than 0.80, except for A. guatemalensis, P. oocarpa, and P. teocote. ^{Baldwin (2009)} mentions that the values of the area under the curve less than 0.7 are considered bad predictive models, values of 0.7 to 0.9 are good models, and values greater than 0.9 are extremely good models. Therefore, 85 % of the distribution models for the pine and fir tree species are considered good (40 %) and extremely good (45 %). Furthermore, the values of the area did not significantly decrease with time and showed a variation coefficient less than 6 % in the majority of cases. The models with low predictive power mainly correspond to those with occurrence localities less than 5, such as A. guatemalensis and P. oocarpa.

Table 3 shows the contribution percentage of the variables or environmental predictors in each species distribution model. The environmental predictors with greater weight were maximum temperature of the warmest month (bio5), minimum temperature of the coldest month (bio6), and pH.

The distribution potential of the 10 pines and two fir trees in the baseline is shown in Figure 2. According to the figure, the suitable area for the distribution of the majority of the species is within the biogeographic provinces of the Trans-Mexican Volcanic Belt and Sierra Madre del Sur; whereas the biogeographic provinces of the Balsas Depression and Pacific Coast are areas where the environmental requirements for the species are not satisfied. Table 4 shows the rate of change of the surface of the suitable area for the 10 pines and two fir trees in future scenarios with regard to the baseline. More than 30 % of the surface area of the state shows ideal conditions for the distribution of P. douglasiana, P. lawsonii, P. montezumae, and P. pseudostrobus. On the other hand, less than 7 % of the surface area of the state is considered ideal for the distribution of A. religiosa and P. rzedowskii. The close outlook predictions show a change in the ideal areas for the presence of the studied species. In this outlook, the distribution area of P. oocarpa will increase 81 %, whereas the suitable area of A. religiosa, P. leiophylla, and P. teocote will decrease by more than 40 %. The predictions for the far outlook also show that there will be a reduction in the ideal areas. The ideal areas of A. religiosa, P. hartwegii, P. leiophylla, P. michoacana, and P. teocote will decrease by more than 60 %. On the other hand, the ideal distribution area for P. oocarpa and P. rzedowskii will increase by 18 and 10 %, respectively.

Figure 2 Potential distribution maps for species of pines and fir in the baseline (1960-2000) of study in Michoacán. 

According to Figure 3, the unsuitable area will increase in the periods and future scenarios with regard to the baseline climate; the increase will cause the decrease of the ideal areas. The decrease of ideal area between the future climate scenarios and the baseline was not significant (P = 0.103), with the exception of the scenario RCP 8.5 and outlook 2075-2099. In general, the species distribution models of the scenario of high emissions (RCP 8.5 W·m^-2) will show a greater reduction of suitable area in comparison with the scenario of low emissions (RCP 4.5 W·m^-2), as indicated by ^{Fernández- Eguiarte et al. (2015b)}. Similar results will be slightly lower in conservative scenarios (^{Gómez-Mendoza & Arriaga, 2007}; ^{Hamann & Wang, 2006}).

Figure 3 Means comparison of the potential distribution area of 12 species of Pinaceae in different periods and scenarios of Michoacán. Different letters in the rows are significantly different in accordance with the Tukey test (P = 0.05).  

The estimation of the decrease of the suitable areas for the presence of A. religiosa is similar to what was found by ^{Sáenz-Romero et al. (2012)}, who reported that the area will be reduced by more than 88 % in periods subsequent to the year 2050. The authors mention that one of the solutions to reduce the impact of climate change on the species distribution is assisted management, locating new localities 300 m above the current altitudinal benchmarks; however, this does not guarantee that the species that coexist (monarch butterfly and edible mushrooms) with A. religiosa will adapt to the new sites.

Due to the increase in temperature and the decrease in precipitation because of climate change, the distribution area of various Pinus species could decrease in areas with altitudes above 3,500 m and it will be necessary to transfer lower altitude species (3,000- 3,500 m) to these altitudinal heights (^{Viveros-Viveros et al., 2009}). ^{Gómez-Mendoza and Arriaga (2007)} found that P. oocarpa is one of the species most vulnerable to climate change in the Mexican republic. According to these authors, the ideal area for P. oocarpa will decrease by 50 %; P. teocote is a species moderately sensitive to climate change, given that its surface area will decrease by 19 %; whereas P. douglasiana is a tolerant species, given that its distribution surface area will only decrease by 1 %. On the other hand, P. leiophylla must be closely studied given that its distribution is not only vulnerable to climate change, but also to some plagues such as Dendroctonus, categorized as a high threat in areas of the Trans-Mexican Volcanic Belt (^{Salinas-Moreno, Ager, Vargas, Hayes, & Zúñiga, 2010}).

The three species that will be most affected by climate change, according to the close and far outlooks, are A. religiosa, P. leiophylla, and P. teocote. These species relate to communities of coniferous and oak forests, whose surfaces could also be affected by the change in soil use due to agricultural, forestry, and livestock activities (Figure 4); therefore, it is estimated that there will be a greater decrease in suitable potential distribution. This is the case with A. guatemalensis, given that the soils where it is commonly found are fertile, an important factor in the establishment of other plantations (^{Andersen, Córdova, Sørensen, & Kollmann, 2006}). Nevertheless, ^{Jantz et al. (2015)} mention that the suitable potential distribution decreased from the year 1500 to 2005, where more than half of the natural vegetation was lost due to such activities. These authors estimate that for future RCP scenarios, more than 70 % of the natural vegetation will be lost due to logging.

Figure 4 Suitable potential distribution of the species that will be most affected by climate change, under the RCP 8.5 2075-2099 scenario. Soil use and vegetation change polygons in accordance with the ^{INEGI (2005)}. 

Conclusions

Of the 27 climate, soil, and topography variables, only seven were utilized (maximum temperature of the warmest month, minimum temperature of the coldest month, annual temperature range, precipitation of the driest quarter, K, pH, and the digital elevation model) with greater significance on the species distribution modeling. The variables generated robust models in time and space, given that 85 % were considered good and extremely good in accordance with the values of the area under the curve. The potential species distribution models in function of the RCP 8.5 scenarios showed the greatest decrease of suitable distribution areas for pines and oaks. The species most affected by climate change will be A. religiosa, P. leiophylla, and P. teocote. Therefore, it is necessary to implement measures to adapt to climate change through assisted management, genetic variation, and plague resistance. The decrease of the ideal areas for the distribution of A. religiosa will be greater than 72 %, being the most affected of the evaluated species; if its population decreases, it will also affect the population of the monarch butterfly and the economy of the neighboring inhabitants who depend on the tourism attracted by the presence of the lepidoptera.