Introduction

The association of alpine tundra with colder climates, and of variation within tundra with various climatic variables, is well established (e.g., Billings, 1988; Walker et al., 1994; Parisod et al., 2010). Many studies have indicated that a warming climate could threaten some alpine tundra (cf. Klanderud and Totland, 2005; Walther et al., 2005; Walker et al., 2006; Lenoir et al., 2008). More specific results on characteristics such as phenology support these concerns (Inouye, 2008; Wipf et al., 2009; Smith et al., 2012). Further, it has long been supposed that forest response to climate warming could figuratively push alpine tundra off the top of mountains (Peters and Darling, 1985), and Diaz and Eischeid (2007) calculated that the climate type associated with alpine tundra will no longer exist in the continental U.S.A. with continued climate change.

The notion of a general loss of alpine tundra has been challenged (Randin et al., 2009a) or at least locally contradicted (Cannone et al., 2008). The geomorphological conditions that create the general variability in microhabitats in alpine tundra are at multiple scales (Lonegran and Del Moral, 1984; Fisk et al., 1998; Liptzin and Seastedt, 2010) and affect the assemblages of alpine tundra in cross-scale interaction. Malanson et al. (2011) reported that the degree of difference in alpine tundra across > 1000 km of the Rocky Mountains could be found within a 4 × 10⁹ m² area in Glacier National Park, Montana (GNP). Mountains have heterogeneous microclimates due to topographic complexity which could limit the impact of regional climate changes (Scherrer and Körner, 2011). The response of alpine tundra to climate change is also variable at other scales (Engler et al., 2011). Our purpose is to examine the impact of geomorphic heterogeneity on variability of alpine tundra plant communities, which is a major context for plant response to climate change. In doing so we assess the ability of simple geomorphic indicators to be proxies for processes of microclimatic modification.

While we agree with Scherrer and Körner's (2011) general conclusions based on their measurement of a wide range of thermal microhabitats within meter-scale distances, we believe that moisture may be as important as temperature generally and more so for many species in alpine tundra. They concluded that “all but the species depending on the very coldest microhabitats will find thermally suitable ‘escape’ habitats within short distances,” but we are concerned that the range of moisture microhabitats may be more spatially variable and that species will not so easily escape drying. Moisture conditions are likely to be an important source of variation in alpine tundra because of the direct need for water (e.g., Harte et al., 1995) and their indirect effect on microbial activity and nutrient dynamics (e.g., Lipson and Monson, 1998). Moisture conditions are likely to be affected by geomorphology at multiple scales (e.g., Litaor et al., 2008). In the Swiss Alps, Vonlanthen et al. (2006) found that the variation in community structure in alpine tundra was best correlated with a gradient of temperature but secondarily with one related to soils, particularly soil tension. In a discussion of topographically-maintained microrefugia, Dobrowski (2011) cited water availability as being as important a limit on species as temperature. Other studies have found effects of topography on diversity (Bruun et al., 2006).

Contributing to moisture variation, snow patterns, especially the timing of snowmelt, have long been known to differentiate alpine tundra (Billings and Bliss, 1959; Walker et al., 1994; Arft et al., 1999; Sandvik et al., 2004; Choler, 2005; Huelber et al., 2006; Bjork and Molau, 2007; Kudo et al., 2010). For example, snow manipulation experiments at Niwot Ridge (Colorado Front Range) changed the nitrogen cycle significantly (Williams et al., 1998). These results have been extended to wider observations where topographically induced variation in snow cover had consequences for species richness (Litaor et al., 2008). Notably, the pattern of snow is affected by topography through its effects on redistribution and melting, and snowbeds are often found in specific relation to topography: on concave leeward slopes and just below leeward ridges (i.e., at cornices) (e.g., Dobrowski, 2011). Geomorphology also channels snowmelt (e.g., Baron et al., 2000; Hood et al., 2003; Choler, 2005), and has a direct effect on tundra where it is associated with a disturbance such as solifluction (Johnson and Billings, 1962; Haugland and Beatty, 2005; Vonlanthen et al., 2008; Randin et al., 2009b).

Here we examine what aspects of environmental variability might allow the development of different kinds of tundra within small areas. We focus on geomorphic variables that modify broader climate factors to produce the microenvironments experienced by plants (which are thus decoupled from regional climates in alpine environments; cf. Pape et al., 2009; Wundram et al., 2010) and thus may serve as a proxy for the more difficult and expensive to measure microclimate and microhydrology variables. Our intent is to provide a context for interpreting change in alpine tundra and a framework for focusing monitoring and mitigation strategies in an era of changing climate. This context is needed because the usefulness of the Global Observation Research Initiative in Alpine Environments (GLORIA; Grabherr et al., 2000) network will depend on observational and theoretical context because prior and ongoing studies have produced variable results (e.g., Harte and Shaw, 1995; Walker et al., 2006; Pauli et al., 2007; Randin et al., 2009a; Abadneh and Woolfenden, 2010). Context is also needed for local understanding and mitigation (Malanson et al., 2006). Alpine tundra is an important amenity resource globally (e.g., Gret-Regamey et al., 2008) as well as a contributor to global biodiversity (Körner, 2003). If we can elucidate the relationship between alpine tundra and local environmental variables that interact with and modify broader climate forcing we will have established a better basis for further monitoring, experiment, and interpretation (cf. Nagy and Grabherr, 2009).

To examine the effects of geomorphology in addition to or in place of climatic variables, and given spatially extensive climate data is lacking in most mountain areas and interpolations are suspect (e.g., Grafius and Malanson, 2009), we will examine the variability of alpine tundra plant community structure in relation to a range of geomorphic indicators with and without broad geographical gradients that may capture coarse scale climatic variability. In this study we make use of existing data that can be an example for other studies. These data include ordinal scale values of species cover on sites for which there are ratio scale transformations, and ordinal scale indicators of abiotic environmental variables. These observations of plant species and their environment are common from older studies (especially in Europe), and potentially can be reused. Such studies can then extend the generality of the ongoing experiments and monitoring that are spatially restricted. In particular, these older studies do not have data on microclimates or the details of water availability. We contend that geomorphic variables mediate macroclimate patterns to produce microclimates, and that within a macroclimate region geomorphic variables can capture much of the variation needed to explain the range of alpine vegetation types.

Study Area and Data Source

Glacier National Park (GNP), Montana, has extensive tundra in geomorphically complex situations (Malanson et al., 2007). Here the tundra is typically floristically diverse. Bamberg and Major (1968) identified 185 species in the area of Siyeh Pass, including the broad expanse of East Flattop Mountain, of which 62 occurred in their sixteen 50 × 20 cm quadrats. Choate (1963) reported 136 species at Logan Pass. She noted that many of the alpine species of GNP are near the southern limits of fairly extensive arctic-alpine ranges. Lesica and McCune (2004) reported that these species were threatened by climate change, while Malanson et al. (2011) found such a broad array of tundra types and habitats to potentially lessen impacts. Lesica (2002) summarized all alpine tundra in GNP into four geomorphic types. One, talus, supports little vegetation due to its instability. What he described as fell fields include relatively rare fell and the broad uplands east of the Continental Divide where solifluction processes often create a stair-step appearance with stony treads alternating with vegetated risers (slightly steeper than the treads but not vertical). Lesica (2002) noted that this vegetation grades into what he called turf (he described these as areas of deeper soil). At the more extreme portion of the gradient is wet meadow. Along this gradient, relatively dry sites are dominated by grasses while sedges dominate the wetter sites. Lastly, Lesica (2002) noted that some types of alpine vegetation are associated with permanent or persistent snow.

Some aspects of the alpine microgeomorphology of GNP have been studied in their own right. At a few isolated convex upland sites in eastern GNP, glaciation did not directly affect mountaintops during the late Pleistocene. Instead, intense periglacial conditions existed that led to the development of turf-banked treads and risers in response to solifluction processes. These treads and risers, described elsewhere (Butler and Malanson, 1989, 1999; Malanson et al., 2002; Walsh et al., 2003; Butler et al., 2004), are the primary environment where alpine tundra is found in GNP. The risers are primarily vegetated by Dryas octopetala and/or Arctostaphylos uva-ursi.

Damm (2001) did the most extensive classification and description of alpine vegetation in GNP. He classified over 500 detailed samples into over 40 associations (a vegetation classification unit based on species composition), some with further divisions, using the European phytosociological system of Braun-Blanquet (1932). He endeavored to use a 16 m² quadrat for his sites, but in places where the tundra was spatially confined he used smaller quadrats. We use these data as representative of the variation of tundra and assume that error would most likely be for rare species that would affect the analyses less than others. Second, we excluded lichens and analyze the patterns with vascular plants and mosses because the two groups respond to environmental gradients in opposite ways (e.g., Hudson and Henry, 2010). In a number of cases, even with good effort, it was not possible to identify all mosses to species, and so only the genus was reported. The rarest species, i.e. those listed in footnotes to the tables, were not included.

Damm (2001) collected basic descriptive information about his sites. These observations include elevation, aspect, and slope—and we use these variables directly, except changing aspect to degrees difference from southwest. His previously unreported field sheets (example included as Appendix 1) include information on topography, geology, soils, surface rocks, and topographic exposure, which we reduce (below). Based on his marking of his field sheets we derived a score of 1–5 for each nominal variable in the category; 1 = no mark; 2 = dotted slash; 3 = solid slash; 4 = solid slash and dotted x; 5 = solid x. Where we had percentage values (e.g., rock cover in size classes) we used them. Recording was sparse for some variables such as soil profile depth and we did not use those. The list of variables for PCA reduction from Damm (2001) is shown in Table 1.

TABLE 1

Variables recorded by Damm (2001) that were reduced using principal components analysis (PCA).

Damm (2001) originally surveyed his sites with a low-resolution GPS, but many sites were resurveyed later with better precision (he had left micro-cairns at his sites and was able to relocate them quite accurately). If Damm (2001) had sampled some types of vegetation all in one area, any test of the importance of site location would be biased. We plotted the UTM coordinates for all sites as a representation of the geographic spread of sampling, and then plotted the coordinates of the plots that Damm (2001) identified as associations and grouped as such in separate tables. In most cases the plots for single associations were scattered across the geographic range of GNP and so it seems clear that site selection was not biased greatly. Given missing data for many cases, however, we finally reduced the number of plots to 396.

For comparison we examined general climate gradients across GNP derived from Daymet data ( http://www.daymet.org). Daymet is a spatial interpolation procedure specifically designed for irregularly spaced stations and complex topography (Hungerford et al., 1989). The available data are from 1980 to 2003 and are at 1 km spatial resolution. We entered the coordinates for alpine locations across the range of GNP to derive the descriptive results (Fig. 1 and Table 2).

TABLE 2

Climate indicators interpolated in DAYMET for points representing the gradients of latitude and longitude across GNP; data are from 1980 to 2003.

FIGURE 1.

The locations of the sites for which Daymet climate data was derived for descriptive context. The latitude and longitude of Logan Pass are 48.6952°N, 113.7180°W.

Methods

Because we had indicators for so many independent variables that relied on Damm's (2001) rankings, we reduced his data on topography, geology, soils, and exposure using principal components analyses (PCA) to derive axes to represent combinations and interactions of these variables for each site. We used PCA separately for the variables in these categories: slope conditions; soils (organic matter and moisture); rock type; exposure; and the cover percentages in the classes of rock size. We used the PCA axes for which the broken-stick eigenvalue exceeded the observed eigenvalue (cf. Frontier, 1976; Jackson, 1993). These PCAs produced two or three axes meeting this criterion that we used subsequently and which we refer to as Slope1–3; Soils1–2; RockT1–2; Expo1–2; and RockS1–3.

We retained Damm's (2001) more direct measures for slope, aspect (which we corrected to difference from 270°), and elevation. We used his geographic coordinates (UTMN, UTME) as variables and derived three additional location-based variables: east or west of the Continental Divide (Eastness: binary, 1,2); Easiness × UTME (EastEast); and Eastness × Aspect (EastAspect). These direct geographic indictors are the factors different from geomorphology and may be related to broader climatic gradients (Table 2 and Fig. 1). The list of independent variables derived by PCA and other direct variables is given in Table 3.

We ordinated the site × species from Damm (2001) using detrended correspondence analysis (DCA) in PC-ORD (McCune and Mefford, 1999); we rescaled the axes to 26 segments and did not downweight rare species. We transformed the Braun-Blanquet cover classes to percent cover using the midpoints of the classes. We also used non-metric multidimensional scaling (NMDS) in preliminary analyses. Ordination provides a mapping of sites in statistical space based on the similarity of their plant community composition. NMDS is a preferred method for many types of exploratory ordination because it makes few assumptions about the form and relations of the data; DCA does make such assumptions and manipulates the data to conform to them, but it often produces results that are more interpretable in terms of relations with environmental gradients. Because NMDS produced a roughly spherical cloud of points, we analyzed the DCA results (with trial rotations we found that a Pearson correlation of an NMDS axis with the primary DCA axis was > .9; we took this result as an indication that the DCA results were not simply an artifact of the method). A combination of ordination methods is a useful quality control approach (e.g., Robbins and Matthews, 2010).

TABLE 3

List of all independent variables.

We next regressed the DCA1 ordination scores on the environmental variables, including the axis scores determined in PCA. We used the “best subset” algorithm in SPSS-19 to select a model, with minimizing the Akaike Information Criterion as the criterion. This linear regression algorithm compares all combinations of independent variables. Because this algorithm does not produce the more widely recognized accumulated R², we also ran a common stepwise multiple regression. Because the direct geographic variables that we assumed correspond to climate gradients are at a different scale (UTMN, UTME, Eastness, EastEast, EastAspect), we ran and interpreted analyses with and without them.

Results

The DCA ordinations produced a spread of points of over 600 and 500 units (100 units is equivalent to one standard deviation of species turnover; Gauch, 1982) on the first two axes (Fig. 2; eigenvalues 0.83 and 0.62, but the latter is not indicative of variance due to rescaling; the coefficients of determination for the correlations between ordination and original n-dimensional distances have cumulative R² of .291 and .385; the overall inertia is 25.59). Given the predominance and better reliability of the first axis (DCA 1), we examined it alone in relation to the environmental variables. The arrangement of plant community types along DCA1 seems to reflect differences in available soil moisture. At one end of the axis are plots classified by Damm (2001) as bearberry communities (northern goldenrod—bearberry association with its bearberry— rough fescue subassociation; based on Solidago multiraiata, Arctostaphylos uva-ursi, and Festuca scabrella) and scree slope communities of Stellaria americana or Saxifraga bronchialis; at the other end are streamside communities of the Saxifraga lyallii variation of the Senecio triangularis—Mimulus lewisii association, drier but deeper soils of beargrass (Xenophyllum tenax) or of Tofieldia glutinosa—Carex lenticularis moss and sedge-dominated plots (Fig. 3).

FIGURE 2.

The pattern of sites on the first two axes of the detrended correspondence analysis (DCA) of the site × species cover data. The plant community types identified are the extremes referred to in the text.

Basic Pearson correlations among all variables indicate a high degree of collinearity in a few cases (only the most relevant are shown here [Table 4]; others are available in Appendix 2). Most notably, given the orientation of GNP and the mountain ranges it encompasses, UTMN and UTME are highly correlated (0.89). Additionally, among the highest correlations are UTME with Rock-Type1 and Exposure1. Due to overall weak support for DCA 2, we do not analyze it further, but our conjecture that the array of community types on Figure 2 could be related to soil moisture is strengthened by the high correlation (0.68) between DCA 2 and Soil1, our best proxy for soil moisture, with the beargrass and moss-sedge communities being at opposite ends of our Soils1 gradient with the bearberry communities in the middle.

The regression analyses produced significant models with good explanatory power (Best Subsets algorithm AIC 3708.7 and Accuracy 70.0%; Stepwise algorithm AIC 3718.8, Accuracy 68.7%, and Adjusted R² .684). The geographic variables UTME and UTMN can account for much of the variance in the models (Table 5), but exposure variables, which are from Damm's estimates of wind and snow duration, are also important.

When we ran these same analyses without the geographic coordinates (noting the collinearity with UTME and UTMN; Table 6), the models are weaker but still significant (Best Subsets algorithm AIC 3822.3 and Accuracy 59.7%; Stepwise algorithm AIC 3826.3, Accuracy 58.9%, and Adjusted R² .575). The exposure variable that captures the degree of protection from wind and the duration of snow now accounts for most of the variance.

Discussion

The reasonable approximation of the models without the geographic variables indicates that the geomorphic indicators used in this analysis are highly related to the variation in alpine tundra plant community composition at the scale of GNP. Without considering the climate gradients within the park, and given that we are not including variables that we expect could be important but for which we have no data (i.e., soil texture and chemistry: e.g., see Bowman and Seastedt, 2001; Loffler and Pape, 2008; or isolation and dispersal: e.g., Marchand et al., 1980; Molau and Larsson, 2000; Vittoz et al., 2009), simple geomorphic indicators account for much of the variation in composition and thus in regional diversity. To the extent to which we do capture soil conditions such as soil moisture, other geomorphic variables explain more variance. Overall, the use of geomorphic variables that are relatively easy to determine in the field proved to be a useful approach in an effort to understand how climate is modified at the plant scale.

FIGURE 3.

Examples of the plant associations found at opposite ends of the primary DCA axis, with characteristic (A) Arctostaphylos uva-ursi (bearberry) and (B) Tofieldia glutinosa—Carex lenticularis (false asphodel—Kellogg sedge).

The topographic variables that are important in the models are those that are associated with available moisture. The PCA axes of exposure represent variation in what Damm (2001) estimated were variations in winds, which affect snow redistribution, and snow cover directly, which includes shade (it is possible that these exposure gradients are related to temperature, since greater surface wind speed should allow greater convectional transfer of heat). The second group of variables includes those of soil characteristics, primarily moisture in Damm's (2001) terms, and rock sizes, that affect drainage and water holding capacity. The third group includes slope and the interaction of slope and elevation, which are key factors in overall soil moisture models for alpine areas (Burns and Tonkin, 1982).

The considerable multicollinearity among the independent variables makes specific attributions of probable cause and effect difficult. For example, the prominence of the geographic variables, especially noting the correlation of UTME and UTMN due to the orientation of the mountains, could be due to their collinearity with other variables (note that the relative importance and R² of the first variable in the models without them is similar to that of UTME alone) but the overall additional explanation could be due to an environmental pattern or artifact (e.g., perhaps the most easterly or westerly areas have unique conditions—but our analysis of the distribution of Damm's (2001) associations would argue against this conclusion), or it could reflect historical processes.

However, topographic variables alone can account for the same amount of variance in the plant community floristic structure as can the addition of geographic coordinates assumed to coincide with major climatic gradients. In the context of Scherrer and Körner's (2011) results, it is notable that the topographic variables that are important either with or without the inclusion of location coordinates are those related to exposure, and defined by Damm (2001) as “wind exposure” and “estimated snow duration.” These are normally the inverse of each other and the PCA axes are both related to the combination of lower winds and longer snow duration in terms of which variables loaded on them. In general, water is considered to be a direct resource for plant growth, while temperature has indirect effects by regulating rates (cf. Austin, 1987).

TABLE 4

Pearson correlation coefficients for selected variables (N = 396; for 2-tailed test: * p < .01; **p < .001).

Topographic factors are often related to snow cover, which may be the proximate cause of resource variation to which plant species are constrained and which differentiates communities. Snow cover may be an environmental factor quite sensitive to climate change. Snow cover would be reduced by warmer temperatures in which a larger proportion of precipitation fell as rain; snow cover would be reduced by higher warm-season temperatures causing faster snowmelt and shorter duration of cover, leading to local drought conditions in late summer; and snow cover can be affected by a change in the temperature and water content of snow as it falls, with warmer, wetter snowfalls being less susceptible to redistribution by wind and so resulting in fewer and/or smaller snow patches, again leading to late summer drought. Those species depend on late summer snowmelt that could be most threatened, complementing those species of the coldest microclimates as identified by Scherrer and Körner (2011). Among those that might be threatened are the plant communities that have the showiest flowers with high appeal to the public, such as the streamside communities of Mimulus lewisii (Fig. 3). Many of these are found at rivulets draining late-lasting snowfields or glaciers that are already waning (Hall and Fagre, 2003; Pederson et al., 2011).

TABLE 5

Linear models developed with the primary DCA axis as the dependent or target variable; direct geographic variables are included.

The importance of topography at multiple scales, and perhaps especially at spatial scales close to the sizes of individual plants, should be of interest in alpine regions more generally, while the particular interaction with snow may be more climate-specific. Glacier National Park experiences snowfall that is often redistributed by wind. Drier snowfalls (e.g., Rocky Mountain National Park, Colorado) may experience redistribution with less local storage, while wetter locales (e.g., Olympic National Park, Washington) may have more stable snowpacks. A focus on the extremes as noted in this study should be considered for monitoring impacts in any case.

TABLE 6

Linear models developed with the primary DCA axis as the dependent or target variable; direct geographic variables are excluded.