This study attempts to model alpine tundra vegetation dynamics in a tundra region in the Qinghai Province of China in response to global warming. We used Raster-based cellular automata and a Geographic Information System to study the spatial and temporal vegetation dynamics. The cellular automata model is implemented with IDRISI's Multi-Criteria Evaluation functionality to simulate the spatial patterns of vegetation change assuming certain scenarios of global mean temperature increase over time. The Vegetation Dynamic Simulation Model calculates a probability surface for each vegetation type, and then combines all vegetation types into a composite map, determined by the maximum likelihood that each vegetation type should distribute to each raster unit. With scenarios of global temperature increase of 1 to 3°C, the vegetation types such as Dry Kobresia Meadow and Dry Potentilla Shrub that are adapted to warm and dry conditions tend to become more dominant in the study area.
Global temperature has increased due to the effects of greenhouse gases emissions, and significant impacts have occurred in the biomes and ecosystem types in the region of the Qinghai-Tibetan plateau (Johnes and Briffa, 1992; Ni, 2000; Leemans, 2004). Based on 38 years (1959–1996) of climate observations and statistical analysis, the annual mean temperature increases during this period ranged from 0.4 to 0.6°C in the area of the Haibei Alpine Tundra Ecosystem Research Station (Li et al., 2004). The equilibrium pattern of plateau biomes in response to climate change was modeled using the BIOME3China model (Ni, 2000). However, depicting the vegetation change in response to global warming is a challenge in terms of using spatial modeling to detect and simulate the vegetation response mechanism to warmer weather. In this paper, we model the alpine tundra vegetation changes by integrating the environmental factors in a raster-based Geographic Information System (GIS) environment. Assuming each vegetation type in the raster cell unit reacts as a homogeneous entity, we conduct a spatial and temporal simulation by integrating cellular automata and Multi-Criteria Evaluation (MCE) provided in the IDRISI software (Eastman, 2003).
Spatial modeling processes are available in current GIS software such as IDRISI, which is capable of dealing with a large set of raster data and manipulating the data via operations in a series of discrete time steps, where single raster cells can be influenced by their neighborhood or other data in an overlay. All map layers are imposed on the same grid system. This type of GIS environment provides a sophisticated tool to help us target the real problem in a complex system (Wolfram, 1984; Coulelis, 1985; Itami, 1994; White and Engelen, 2000; Giles, 2002). In our study, we use GIS analysis, linear regression, MCE, cellular automata (CA), and a raster image calculator to build a unique Vegetation Dynamic Simulation Model (VDSM). Global warming scenarios are interpreted as inputs of the spatial parameters. Large processing tasks are completed by the computer system.
The predicted outcome of this study is that individual vegetation types will respond to a global mean temperature increase (GMTI) in 2100 of 1 or 3°C by either expanding or shrinking their range because of plant species' suitability to the warmer and drier climate conditions. This corresponds to 0.1 or 0.3°C per decade, respectively (Leemans, 2004).
Our study site is located near the Haibei Alpine Tundra Ecosystem Research Station, Qinghai Province, China (37°29′–37°45′N, 101°12′–101°33′E) (Zhang and Zhou, 1992). The elevation in the study area varies from 3000 to 4500 m a.s.l. The model uses the following data:
A 90 m × 90 m resolution DEM;
Temperature derived from the DEM using an empirical linear regression model (Zhang, 2005);
Land surface parameters (aspect, slope, stream channel density) derived from the DEM using GIS analysis tool provided in the IDRISI software (Eastman, 2003);
Constraints are the limited areas that are not considered to be natural vegetation, such as water bodies, glaciers, gravel slopes, and artificial grasslands. They represent areas where the natural vegetation cannot grow or are otherwise constrained. These areas are excluded from the vegetation map. A Boolean image is created to display inclusion and exclusion of the constraint conditions.
The factors (F1, F2, …Fn) used in the MCE are selected based on the most important variables that determine the output y in Equation (1). In our study, we use MCE in step 1 (Fig. 2) to determine a normalized surface temperature, which calibrates the temperature by these spatial parameters: aspect, suitable surface temperature, and distance to the nearest stream channel.
The temperature has a certain spatial pattern in the study area: (1) It increases from the north (a = 0°) to the southeast (b = 145°), with the highest values from the southeast to the southwest (c = 275°), then decreasing to the north (d = 360°). The change pattern can be described as a sigmoidal fuzzy function type, in a symmetric shape with specific values (a, b, c, d) in Table 1. (2) The distance to a stream affects the temperature. If the area is within 10 m of a stream or water body, its temperature is closer to the stream or water temperature; beyond 600 m from a stream or water body, the temperature is minimally affected by the nearest water bodies. The values are defined as a = 10 m, b = 600 m, with a sigmoidal fuzzy function type and monotonically deceasing shape. (3) In the lower valley of our study area, the monthly mean temperature in July is 10.1°C (Li et al., 2004). The temperature decreases with increasing elevation. We define the temperature less than a = 0°C as unsuitable for alpine plant growth. The temperature from 0°C to b = 5°C is defined as less suitable for alpine plant growth; the temperature from 5°C to c = 13°C is defined as most suitable. The temperature from 15.5°C to d = 18°C reduces the suitability for alpine plant growth, representing the dry, south-facing areas (Table 1). In the Analytical Hierarchy Process (AHP), the most important temperature is in the range of b = 5°C to c = 13°C, which present the sigmoidal fuzzy function type and symmetrical shape.
Normalized surface temperature's fuzzy membership function types and shapes, factor's suitability values (a, b, c, and d), and AHP weights.
In step 2 (Fig. 2), we use the MCE method and create a suitability map for each vegetation type using the aspect, slope, and the distance to stream channels. The suitability values (a, b, c, d) in Table 2 are defined based on the vegetation distribution in the study area (Zhang and Zhou, 1992). The fuzzy membership function type and shape and AHP weights for each vegetation are calculated and reported in Table 2 (Eastman, 2003).
Ten vegetation types' fuzzy membership function types and shapes, factor's suitability values (a, b, c, and d), and AHP weights.
The macro model is created to simulate changes in each vegetation type through time and space, integrating operations such as overlay, scalar, fuzzy module, and cellatom (Fig. 3). These operations are available in the IDRISI software package (Eastman, 2003) and have to be built in Macro Modeler with an initial scalar value (0.0–1.0) and fuzzy set values (0–255).
The GMTI scenarios are implemented in the Macro Modeler by adjusting the scalar operation to increase the temperature by 0.1°C in a discrete time period. The same logic is applicable to 0.3°C in a discrete time period (Leemans, 2004). Running the simulation for 10 iterations, the effects of increasing temperature on each vegetation type are accumulated in the output image.
The cellular automata module is implemented in the Macro Modeler (Fig. 3) and used to form a uniform raster image to represent global warming effects in a spatial context, which operates over discrete time steps (Coulelis, 1985; Giles, 2002; White and Engelen, 2000). The change in each cell depends on the parameters or requirements set by the user and the surrounding neighbors (Wolfram, 1984; Itami, 1994; Ruxton and Saravia, 1998). This project used the CELLATOM module with a 3*3 filter and reclassifies an output cell if at least 3 neighbors contain non-null values. We define the 10 iterations using the DynaLink module (Eastman, 2003). Within each iteration, each vegetation map is dynamically updated by running Cellatom, then it is overlaid with the vegetation suitability map altered by a GMTI of 0.1 or 0.3°C. Thus, after 10 iterations, a final suitable vegetation map is produced. The same dynamic processing is repeated for each vegetation type resulting in a total of 10 vegetation suitability maps in response to warmer weather.
Composite Final Vegetation Map
Our objective is to create a composite vegetation map for each global warming scenario, GMTI of 1 or 3°C over time (Fig. 2). All of the 10 vegetation suitability maps with a GMTI of 1 or 3°C are combined in order to produce a composite map using the image calculator module in IDRISI (Eastman, 2003). The highest suitability among the vegetation types is selected to represent the successful vegetation type in that cell.
Equation (2) above creates a map where the value at each pixel corresponds to the vegetation type with the highest probability of thriving, where Veg1, Veg2 … Veg10 correspond to each of the 10 vegetation types. The 10 resulting maps of dominant vegetation types are overlaid to produce a composite distribution of the most probable vegetation types.
Normalized Temperature Spatial Distribution
Temperature changes across the study area are not only due to elevation, but also due to aspect and distance from the nearest stream channel. The linear regression model provided a temperature spatial distribution based on elevation alone, which is our primary step. Furthermore, the normalized temperature surface created by the MCE is highly representative of the potential temperature distribution in a normalized fuzzy format (Fig. 4).
Vegetation Change Comparison
We calculate the percent area change for each vegetation type (Fig. 5) before we compose the final vegetation map. By increasing the global mean temperature (0.1 or 0.3°C per decade), the percent area change of each vegetation type indicates the potential expansion from their original ranges. For example, without competing with other types of vegetation, both Dry Potentilla Shrub and Dry Kobresia Meadow could expand their vegetation area by 23%. The difference between them is that Dry Potentilla Shrub responds more positively to a GMTI of 1°C, and Dry Kobresia Meadow responds more positively to a GMTI of 3°C. Similar phenomena are also observed in other vegetation types.
After we compose the final vegetation map, the highest suitability among the vegetation types is selected to represent the successful vegetation type in every cell. For instance, Dry Potentilla and Dry Kobresia expand into areas previously occupied by wet types of vegetation (Fig. 6). The Riverside Blysmus meadow, which requires moist conditions, disappears completely with a 3°C temperature increase. In general, the dry vegetation types demonstrate significant expansion from their original ranges and tend to become more dominant in the study area.
Vegetation Dynamics Over Time
The time dimension in the CA module of IDRISI is interpreted as a discrete time step, which corresponds to a temperature-time dimension. The CA module is integrated into the macro modeler, and uses the DynaLink module (Eastman, 2003) to simulate vegetation change within each iteration. The GMTI is defined as a spatial parameter in the deterministic model with a temperature increase of 0.1 or 0.3°C per decade.
Vegetation Dynamic Simulation Model (vdsm)
The VDSM is an example showing that spatial modeling can solve complex ecosystem problems in terms of having the capability to simulate spatial and temporal vegetation dynamics. In this case study, we model the spatial distribution of temperature and create GMTI scenarios as a spatial grid image. The vegetation dynamics are simulated in discrete time by applying CA in a Macro Modeler. In future studies, this model will be capable of modeling the water-time dimension that makes the simulation more adaptable to global warming research. The VDSM could be potentially incorporated with a normal climate change model to assist in a long-term ecosystem simulation. Alternatively, the VDSM is capable of linking with a stochastic model of temperature change, in which we might be able to forecast an ecosystem disaster.
Temperature distribution is correlated with and controlled primarily by elevation. Numerous spatial interpretation methods have been applied to estimate the spatial distribution of temperature (Li, 2005). The interpolation results do not always agree with the actual sample points, including using geo-statistical methods, and spatio-temporal spline. These methods are highly dependent on the distance to the sample points, and the surface equation. In our study, our first step is to create the primary temperature surface based a linear relation with elevation. The objective is to obtain a more accurate temperature map in terms of aspect, suitable temperature, and distance to the stream. We use the Multi-Criteria Evaluation with Weighted Linear Combination (MCE_WLC) to calibrate the spatial temperature distribution. The fuzzy memberships between the temperature and each factor (aspect, suitable temperature, distance to stream) are based on previous research works (Zhang and Zhou, 1992; Zhang and Welker, 1996; Zhang, 2005). The output, the normalized temperature surface is set into fuzzy format (0–255). Since the temperature is major factor on determining vegetation composition, structure, and distribution, the normalized temperature surface plays an important role when we simulate vegetation dynamics in spatial and temporal dimensions.
Error and Fuzzy Sets
Data from multiple sources are unlikely to be the same resolution, scale, or nominal context. The errors occurred during data conversion, and propagated through the simulation. This model did not take competition among plant species and grazing effects into account (Zhang, 1990; Wu et al., 1996). Clearly, some species are more aggressive and dominant than others, and this would affect their probability of survival when compared to other less aggressive plant species. The lack of a grazing factor resulted in some uncertainty in the modeled vegetation changes.
Fuzzy membership function types and shapes for the factors are based mainly on previous research. We had difficulty validating the accuracy of values of a, b, c, d on Tables 1 and 2. But, the weights of AHP illustrate the importance of each factor in terms of determining the vegetation distribution. For example, the Wet Potentilla Shrub is distributed on north-facing slopes, where the north aspect is more important than the slope and the distance to a stream (e.g., aspect weight 0.5842 > slope weight 0.2318 > distance to stream weight 0.184). Similarly, the Wet Willow Shrub, Caragana Shrub, Wet Kobresia Meadow, Dry Kobresia Meadow, and Carex Meadow are distributed and correlated with the different aspects of the terrain. However, the Dry Potentilla Shrub and Typical Kobresia Meadow are distributed on open valley, and the aspect and slope become less important to their distributions. Also, the Riverside Blysmus Meadow and Riverside Kobresia Meadow are mainly distributed between 3200 and 3300 m within a buffer area 600 m from the stream or water bodies. As a result, the spatial factors' fuzzy sets agree strongly with the vegetation distributions.
Vegetation Dynamic Simulation Model
The VDSM starts from thinking about how to solve a complex vegetation dynamics problem using CA. The VDSM is built by combining the MCE, Macro Modeler, CA, image calculator, Scalar, and Fuzzy functions in IDRISI. Making a clear objective helps us to look into many available modelers and functions in order to solve the problem within the IDRISI software environment. For instance, using Decision Wizard, we create an objective and define a set of constraints, which eliminates the areas that are not natural vegetation. The spatial patterns of the factors (temperature, aspect, slope, and distance to stream) are created as continuous surfaces using Fuzzy functions. The transition rules in Macro Modeler are defined as the maximum potential suitable vegetation in each cell as well as over the study area. The frame work of VDSM (Fig. 2) is a summary of the model structure and functionalities. VDSM is flexible enough to be integrated with other sub-models that are available from GIScience technology. Figure 3 provides an example of Macro Modeler incorporating Cellular Automata for the modeling of Caragana Shrub. In the case of a GMTI of 0.1°C, the suitability map for Caragana is weighted at 0.4 and the normalized temperature map is weighted at 0.6. These two factors are combined using the overlay module to produce a map of suitability for Caragana Shrub with an incremental temperature increase of 0.1°C. Ten iterations on one vegetation layer are simulated and updated after each iteration using the Cellatom and Dynalink modules. At the end of the iterations, we obtained an accumulated effect of GMTI of 1°C on the vegetation layer.
It would also be possible to incorporate the water layer with the vegetation layer and temperature layer, which can be linked by its weighted value in the VDSM. Thus, the VDSM not only provides a discrete time representation, but also demonstrates how we could develop our model for use with more complex scenarios (Chapin et al., 2006; McGuire et al., 2006).
Compositing a final vegetation map demonstrates the power of GIS analysis in IDRISI. The image calculator module in IDRISI successfully carries out interpolation of the image calculation presented in Equation (2) at the grid cell level. The VDSM illustrates how we could study vegetation dynamics and model many other spatio-temporal phenomena.
The VDSM integrates the suitability maps created from MCE, Macro-Modeler, CA, and spatial environmental factors. And the temperature-time dimension model is incorporated into the VDSM, which makes the temperature a spatial parameter that affects the vegetation dynamics over a discrete time step. The simulating processes conducted by Macao Modeler generate the temperature increase of 0.1 to 0.3 °C per decade, which represents the influences of the different global warming scenarios. The results from Figures 5 and 6 demonstrate that global temperature increase reduces moisture availability (Zhang and Welker, 1996) such that dry vegetation can invade areas previously occupied by vegetation adapted to moist conditions. The structure of the model is generally applicable to other situations, but the particular factors and constraints used in this model are unique to the Haibei alpine tundra ecosystem.
Global warming has strong effects on the alpine ecosystems in terms of altering the biomes and ecosystem biodiversity (Zhang, 1993; Ni, 2000; Song et al., 2005). The alpine ecosystem in the region of the Qinghai-Tibetan plateau is sensitive and vulnerable to the changing climate (Zhang and Welker, 1996). The VDSM illustrates that altering global mean temperature changes the alpine vegetation dynamics. With the future integration of water condition (Hodkinson et al., 1999) and disturbance regimes (Zhang, 1990; Zhang and Liu, 2003; Chapin et al., 2006) into the VDSM, the simulation could model more detailed mechanisms and complex feedbacks (McGuire, 2006) of the alpine tundra ecosystem to the changing climate.
We thank Dr. Dragicevic, Verda Kocabas, and the Simon Fraser University Spatial Information Systems (SIS) lab for their support of this research project. We thank the two reviewers who gave us the critical suggestions to improve our revised manuscripts. We thank the Northwest Plateau Institute of Biology, the Chinese Academy of Sciences, for academic support and encouragement of this project.