Arthropod communities along an elevation gradient in Denali National Park and Preserve, Alaska: Rapidly shrinking tundra hosts a unique assemblage of specialists

ABSTRACT

Arthropods at high latitudes and elevations are likely to be vulnerable to effects from climate change such as increased temperatures and shifting vegetation boundaries. Though range shifts northwards and upslope have been reported for many arthropod taxa in temperate latitudes, baseline data needed to track such changes are scarce at northern latitudes. We investigated the influence of climate and vegetation cover on the abundance, diversity, and species composition of pollinators and epigeic arthropods along an elevation gradient in Denali National Park and Preserve, Alaska. We compared arthropods across three habitat types: low-elevation forest, mid-elevation shrubs, and high-elevation tundra along five replicate transects over three years. We collected 35,473 arthropods representing 510 species. Arthropod communities differed distinctly across the three habitat types, with tundra having the highest number of strong indicators and unique species. Elevation, air temperature, and vegetation structure were strong drivers for the ordination of sites. As treeline and shrubline shift upslope with climate change, we predict that distributions of some arthropods will shift to track these habitat boundaries and that tundra-associated arthropods will be most vulnerable as their habitat shrinks. Long-term monitoring of arthropods along elevation gradients at northern latitudes is needed to detect such declines.

KEYWORDS:

• Pollinators
• insects
• elevation
• tundra
• climate change

Introduction

Dramatic declines in arthropod biomass and species richness have been documented around the globe (Dirzo et al. 2014; Hallman et al. 2017; Janzen and Hallwachs 2019). Though causes for individual declines can be challenging to determine, there is general consensus that a multitude of anthropogenic threats are causing significant changes to the abundance, diversity, range limits, and other attributes of arthropods (Sánchez-Bayo and Wyckhuys 2019; Wagner et al. 2021).

Among these threats, climate change stands out as a pervasive and accelerating global phenomenon, affecting arthropods in both human-dominated landscapes and remote wilderness regions (Soroye, Newbold, and Kerr 2020; Halsch et al. 2021; Raven and Wagner 2021). Arthropods in high-latitude and high-elevation environments are especially vulnerable to warming temperatures and associated changes such as advanced snowmelt and changes in precipitation (Parmesan 2006; Inouye 2020). High latitudes have a warming rate that far exceeds the global average (Serreze et al. 2000; Ballinger et al. 2020), and arctic and alpine tundra biomes are among the most sensitive to climate change as trees and shrubs encroach from lower latitudes and elevations (Myers-Smith et al. 2011; Wang et al. 2012; Mekonnen et al. 2021).

Range shifts northwards and upslope in response to climate warming have been reported for many arthropod taxa, including spiders, dragonflies, butterflies, bees, and ground beetles (Hickling et al. 2006; Wilson et al. 2007; Chen et al. 2009; Pyke et al. 2016). However, a landmark study by Kerr et al. (2015) showed that over the last century bumble bee species across North America and Europe have suffered range compression due to contraction at their southern limits without expansion at their northern limits. Additionally, northern species did not retreat to higher elevations as was the trend for southern species. Clearly, predicting how climate change will affect ranges of northern and alpine arthropod taxa requires consideration of more than thermal limits, due to idiosyncratic differences in life history, dispersal ability, extinction risk, and other traits (Hodkinson 2005; Parmesan 2006; Franzén and Öckinger 2012).

Investigating patterns of arthropod distribution along elevation gradients can provide insights on relationships between climate and species distributions, which in turn can inform predictions about future species responses to climate change. Elevation gradients in high-latitude mountains are spatially compressed so that vegetation transitions from boreal forest to alpine tundra may occur within several hundred meters of elevation change (Roland and Schmidt 2015). In Alaska, there is mounting evidence that elevation range limits of some animals are tracking climate-driven upward shifts in treeline and shrubline (Lloyd and Fastie 2003; Dial et al. 2007; Mizel et al. 2016). As yet, we do not have a clear understanding of how arthropod communities differ across habitats along these compressed elevation gradients or how arthropods will respond to changing vegetation boundaries, but shifts in community composition and dynamics could have significant cascading effects in both boreal and tundra ecosystems.

Arthropods are the most abundant and diverse component of the montane fauna (J. S. Edwards 1987) and perform key functions such as pollination, predation, and decomposition that are integral to all montane environments supporting life (Hodkinson 2005). Arthropods also possess numerous traits that make them potentially powerful bioindicators on montane gradients: many high-elevation species have narrow thermal niches, and their metabolism has adapted to relatively low temperatures (Mani 1968); short life cycles and high dispersal ability allow some arthropods to closely track their optimal climatic niches (Danks 2004); and they are often sensitive to fine-scale changes in vegetation cover (Hansen et al. 2016). Though measuring arthropod abundance and diversity along elevation gradients has been a topic of extensive research (e.g., Kearns 1992; Bowden and Buddle 2010a; Viterbi 2013), finding patterns that hold true across multiple taxa, sampling scales, and geographic regions has proven to be somewhat elusive (McCoy 1990). Thus, there remains an urgent need for baseline and long-term elevation studies on a range of arthropod taxa in far northern mountains, where change is occurring at a rapid and accelerating rate.

Our study investigated the abundance, diversity, and species composition of pollinators and epigeic arthropods along an elevation gradient in Denali National Park and Preserve (DNPP) at the southern boundary of interior Alaska (63° N latitude). Pollinators were represented by bees (Hymenoptera: Anthophila) and syrphid flies (Diptera: Syrphidae). Syrphid flies are twice as diverse as bees in Alaska and their free-living larvae include predators, herbivores, and detritivores (Skevington et al. 2019). The epigeic taxa we targeted were beetles (Coleoptera), spiders (Araneae), and grasshoppers (Orthoptera: Acrididae). These taxa represent a diverse array of feeding habits, dispersal capabilities, and sensitivity to abiotic factors, including microclimate. Both pollinators and epigeic arthropods have been used worldwide as bioindicators in ecological studies assessing habitat partitioning along ecological gradients (Sommaggio 1999; Rykken, Moldenke, and Olson 2007; Bowden and Buddle 2010a; Høye et al. 2018; Kazenel et al. 2020).

We sampled arthropods in three topographically defined vegetation communities (referred to as habitat types in this study) along a montane elevation gradient in DNPP to represent successively higher elevation zones: low-elevation boreal forest, mid-elevation tall shrubs, and high-elevation alpine tundra. Specifically, our objectives were to (1) determine whether arthropod abundance, diversity, and community composition differed between habitat types along the elevation gradient; (2) identify strong indicator species associated with each habitat type; (3) examine relationships between observed arthropod distributions and climatic and measured vegetation variables; and (4) establish a robust sampling design, baseline database, and reference voucher collection of specimens that can be used for future long-term monitoring of arthropods in DNPP.

Materials and methods

Study area

Denali National Park and Preserve encompasses 2,452,400 ha in south-central interior Alaska. DNPP straddles the central portion of the Alaska Range, which includes the tallest mountain in North America, Denali (6,190 m). Our study area was located in foothill ranges on the north side of the Alaska Range, typified by mountains with steep alpine talus slopes and lower hills flanked by Pleistocene glacial deposits, underlain by discontinuous or sporadic permafrost (Clark and Duffy 2006). The park road runs east–west for 148 km through this section, allowing easy access to both lowlands and uplands.

The climate north of the Alaska Range is classified as subarctic continental (Clark and Duffy 2006), with long, cold, dry winters and short, warm summers. The snow-free season at high elevation lasts from June through August; these are the warmest and wettest months in DNPP. During the period 1981–2010, the mean summer air temperature (June–August) near the park entrance (63.6500 N, −148.8000 W; 625 m elevation) was 11.7°C, and the precipitation total was 205.2 mm (Sousanes and Hill 2020). Mean summer air temperature in DNPP ranged between 0.2°C and 1.4°C above this long-term average during the years we collected arthropods and/or temperature measurements (2016–2019; Sousanes and Hill 2016, 2017, 2018, 2019).

Forested areas are composed primarily of white spruce (Picea glauca), which peaks in abundance at 450 to 500 m and declines with elevation to an upper limit of about 1,100 m (Roland, Schmidt, and Nicklen 2013). Black spruce (Picea mariana) occurs in lowland wet, acidic sites, and a few broadleaf species, quaking aspen (Populus tremuloides), Alaska birch (Betula neoalaskana), and balsam poplar (Populus balsamifera), are associated primarily with disturbed areas and south-facing slopes (Roland, Schmidt, and Nicklen 2013). Tall shrubs, primarily willow (Salix spp.), birch (Betula spp.), and alder (Alnus viridis), peak in the subalpine zone between 800 and 950 m (Roland and Schmidt 2015). True alpine tundra vegetation, including dwarf shrubs, graminoids, and forbs, extends from around 950 m to the perennial snow line (Roland and Schmidt 2015). The alpine zone also encompasses rocky fellfields and barren scree and talus slopes.

Soil moisture and organic depth generally decrease with increasing elevation, and pH is lower in the lower elevation organic soils and higher in the more mineral alpine soils (Roland and Schmidt 2015).

Site selection

We established five sampling transects along the park road in the northeastern portion of DNPP (Figure 1). Each transect consisted of a low-elevation forest site, a mid-elevation tall shrub site, and a high-elevation alpine tundra site (Figure 2). Thus, we sampled five replicate sites for each of three habitats, totaling fifteen sites in all. Low-elevation forest sites were dominated by Picea glauca, either monotonic stands or mixed with P. mariana and Populus tremuloides. The forest understory was generally open, with herbaceous plants such as Linnaeus borealis, Mertensia paniculata, Rhododendron groenlandicum, Rosa acicularis, Vaccinium spp., and Valeriana capitata and a thick ground cover of mosses. Mid-elevation shrub sites were dominated by tall shrubs >1 m in height, predominantly Alnus, Betula, and Salix species. The understory in shrub sites included Aconitum delphinifolium, Mertensia paniculata, Potentilla fruticosa, Rhododendron groenlandicum, and Vaccinium vitis-idaea. High-elevation tundra sites were dominated by low-stature herbs (e.g., Aconitum delphinifolium, Arnica spp., Gentiana spp., Myosotis alpestris, Pedicularis spp., Polygonum bistorta, Sedum rosea, Silene acaulis), graminoids (e.g., Carex and Poa spp.), dwarf shrubs (e.g., Cassiope tetragona, Diapensia lapponica, Dryas octopetala, Rhododendron spp., Salix spp., Vaccinium spp.), mosses, and lichens (Hultén 1968).

Figure 1. The study area in Denali National Park and Preserve, Alaska. We sampled along five elevational transects (1–5), each containing a forest (F), shrub (S), and alpine tundra (T) site.

Figure 2. An example of the stratified vegetation gradient ascending a mountain side in Denali National Park and Preserve shows the locations of forest (F), tall shrubs (S), and alpine tundra (T).

Transect locations were nonrandomly selected to allow easy access on foot (< 4 km) to all three habitat types from the road. Sites within a transect were spaced 1,300 to 8,400 m apart and at least 100 m from the road to minimize the influence of dust and edge effects. Elevation of sites ranged between 595 and 1,252 m (Table 1).

Table 1. Location of study sites in Denali National Park and Preserve.

Site	Habitat type	Locality	Latitude (N)	Longitude (W)	Elevation (m)
1F	Forest	Mt. Healy	63.73478	148.92911	595
1S	Shrub	Mt. Healy	63.74016	148.95346	927
1T	Tundra	Mt. Healy	63.75095	148.96071	1,183
2F	Forest	Sanctuary River	63.71967	149.48056	772
2S	Shrub	Primrose Ridge	63.72686	149.40852	883
2T	Tundra	Primrose Ridge	63.74163	149.35381	1,252
3F	Forest	Teklanika River	63.64261	149.57181	833
3S	Shrub	Igloo Canyon	63.59455	149.60834	970
3T	Tundra	Igloo Canyon	63.57143	149.64784	1,200
4F	Forest	Toklat Road Camp	63.52373	150.04829	960
4S	Shrub	Caribou Creek	63.50232	150.06764	1,047
4T	Tundra	Caribou Creek	63.49723	150.10019	1,179
5F	Forest	Wonder Lake	63.48534	150.86235	674
5S	Shrub	Wickersham Dome	63.54800	150.97749	811
5T	Tundra	Wickersham Dome	63.55630	150.97755	977

Note. Geocoordinates in datum WGS84.

Arthropod sampling

We selected two pollinator taxa for our study: bees (Hymenoptera: Anthophila) and flower flies (Diptera: Syrphidae). Though these two groups comprise only a subset of the pollinator fauna in DNPP, they represent two of the most conspicuous, diverse, and effective pollinator taxa. We selected three epigeic arthropod taxa: spiders (Araneae), beetles (Coleoptera), and grasshoppers (Orthoptera: Acrididae). Beetles and spiders are among the most abundant and diverse taxa in high-latitude ecosystems (J. S. Edwards 1987; Danks 1993). Grasshoppers, though far less diverse than beetles or spiders, were included to represent herbivores because high-latitude beetle communities are dominated by carnivores (Ernst and Buddle 2015). Grasshoppers are also among the largest-bodied arthropods in DNPP and relatively abundant, thus making them a likely valuable food source for insectivorous vertebrates. All selected arthropod taxa are relatively well known taxonomically and can be efficiently sampled with standardized protocols.

We collected arthropods in 2016–2018 using the sampling methods described below. Arthropods were sampled at two-week intervals across three seasons—2 July to 18 August 2016 (three samples); 13 June to 11 August 2017 (four samples); 5 June to 3 August 2018 (four samples)—for a total of eleven sampling periods over three years. In 2016, we sampled only transects 1 to 4. In 2017 and 2018, we established an additional transect (5) at the west end of the park road (Figure 1).

To sample epigeic taxa, we used a combination of pitfall trapping and litter sampling. At each site, four pairs of pitfall traps were deployed in the corners of a 25 m × 25 m grid. Paired traps were spaced 1 m apart. Each trap consisted of two nested yellow plastic cups, 10 cm in diameter. The cups were buried with the mouth of the uppermost cup flush with the ground surface. We filled each trap one-third full of propylene glycol (Peak Sierra antifreeze) and added ~1 mL of denatonium benzoate as a bittering agent to deter curious vertebrates. Uncovered yellow cups may collect a more diverse sample by acting as both pitfall and pan trap (Høye et al. 2018), but open cups can fill with rain and debris, so one trap from each pair was randomly selected and covered with an aluminum roof. The traps were emptied at the end of each sampling period by straining contents through a fine mesh cloth, which was then inserted in a Whirl-Pak plastic bag for transport to the laboratory. The propylene glycol was replaced for each sampling period and waste propylene glycol was removed from the field.

To sample litter arthropods, we collected approximately 0.01 m³ of leaf litter haphazardly from each site. The litter was sifted to remove large particles and to increase the specimen density of each sample (Guénard and Lucky 2011). The samples were stored in a canvas bag for up to three days before being transferred to a Winkler or Berlese funnel for extraction. Extraction proceeded for fourteen days as suggested by Owens and Carlton (2015).

We used bee bowls (pan traps) and blue vane traps to collect pollinators. Bee bowls were plastic cups (6 cm diameter at the mouth and 10 cm deep) painted white, yellow, and blue. Each cup was raised approximately 20 cm off the ground by a wire support stake. Bee bowls were set out in transects of six or fifteen cups, spaced 5 m apart with alternating colors, and centered in the middle of the site. Each cup was filled approximately one-third full with a 1:1 solution of water and propylene glycol (Peak Sierra antifreeze). At each collection we pooled samples from all traps by straining trap contents through a tight-mesh kitchen strainer and transferring the specimens into a Whirl-Pak plastic bag with 95 percent EtOH.

Blue vane traps are known to attract a wide variety of bees as well as other pollinators (Stephen and Rao 2005; Hall and Reboud 2019). Each vane trap (www.bluevanetraps.com) consisted of a 1.9-L yellow plastic jar with a blue 15-cm-diameter screwcap funnel and two blue 24 cm (length) × 13 cm (width) plastic cross vanes inserted in the top of the funnel. Vane traps were attached to plastic posts approximately 1 m above the ground. One or two vane traps were located at a site. Vanes were separated by 50 m when two were deployed simultaneously in one site and 10 m from the nearest bee bowl. Each vane trap was filled with 3 cm of a 1:1 solution of water and propylene glycol. At each collection, samples were extricated and preserved in the same manner as for bee bowls.

Bee bowls and vane traps were deployed in various combinations (zero to fifteen bee bowls; zero to two vane traps) for different sampling periods over the three years of the study (due to field logistics and pollinator phenology; e.g., minimizing bumble bee queen captures early in the season). However, for each sampling period, exactly the same combination and number of traps were used at every site. Ultimately, one vane trap and six bee bowls were found to collect adequate samples at each site throughout the season and this sampling design is recommended for future work.

All specimens were sorted to morphospecies under a dissecting microscope. Mature individuals were identified to species when possible and left as morphospecies when not. Identifications were determined by the authors using published and online keys (bees: Sheffield et al. 2011; Williams et al. 2014; Ascher and Pickering 2019; beetles: Hatch 1957; Campbell 1968; Lindroth 1969; El-Moursy 1970; O’Brien 1970; Lawrence 1971; Wilcox 1972; Campbell 1973; Bright 1976; Campbell 1979, 1982a, 1982b, 1983, 1984; Anderson and Peck 1985; Gordon 1985; Johnson 1986; Smetana 1988; Campbell 1991; Johnson 1991; Thomas 2000; Arnett et al. 2002; Bright and Bouchard 2008; Pelletier and Hébert 2014; Bousquet et al. 2017; Pelletier and Hébert 2019; grasshoppers: Catling 2008; spiders: Kaston 1946; Chamberlin 1948; Levi 1957; Chamberlin and Gertsch 1958; Levi 1962, 1963; Leech 1972; Opell and Beatty 1976; Dondale and Redner 1978, 1982, 1990; Platnick and Dondale 1992; Dondale et al. 2003; Ubick, Dupérré, and Roth 2005; Whitman-Zai et al. 2015; Marusik, Omelko, and Ponomarev 2017) or by taxonomic experts (see acknowledgments). All specimens were vouchered into the University of Alaska Museum Insect Collection, with the exception of a synoptic pollinator collection retained by DNPP. The specimen data are available publicly through the online database Arctos, which serves data to Global Biodiversity Information Facility and other online aggregators.

Environmental variables

To describe the ground cover within sites, we conducted vegetation surveys in the summer of 2018 at five random locations within each site. A 1-m² quadrat frame was laid on top of the understory vegetation and then percentage cover of plant growth forms (listed in Table 2) was estimated within the frame. These variables were then reduced by principal component analysis to create new orthogonal variables (Table 2). We retained those components that, combined, explained at least 80 percent of the observed variance for downstream analysis.

Table 2. Principal component loadings of ground cover types recorded in fifteen sites across forest, shrub, and tundra habitat types in Denali National Park and Preserve.

Ground cover types	PC1	PC2	PC3
Evergreen trees	3.54	−1.01	−0.04
Deciduous tree	0.25	−0.91	−0.47
Evergreen tall shrubs (>1 m)	0.41	−0.04	0.03
Deciduous tall shrubs (>1 m)	−0.61	−0.45	−0.53
Evergreen low shrubs (20–100 cm)	0.10	−0.01	0.01
Deciduous low shrubs (20–100 cm)	−1.21	5.16	−2.31
Evergreen dwarf shrubs (<20 cm)	0.37	−0.18	0.11
Deciduous dwarf shrubs (<20 cm)	−1.75	0.70	−1.12
Evergreen prostrate shrubs	−0.06	0.35	0.13
Deciduous prostrate shrubs	−0.67	−1.48	−0.10
Forbs	0.15	−0.21	−0.16
Grass	−0.40	−0.23	−0.35
Sedge	0.38	−0.02	0.14
Moss	9.33	1.63	0.60
Lichen	−2.27	2.26	4.73
Bare soil	−0.04	−0.03	0.01
Rock	−0.23	−0.16	−0.24

Note. The first three principal components explained 80 percent of the variation among sites collectively. The most influential eigenvalues are in bold.

We measured air temperature hourly at the center of each site with HOBO 8 K Pendant Data Loggers. Data loggers were deployed 10 June to 15 August in 2018 and 2019. The 2019 temperature data were not concurrent with arthropod sampling but were included to account for interannual variability. Each air temperature logger was placed inside a solar radiation shield (Onset RS1) that was fastened to a metal rod at 1.5 m above the ground. To compare air temperature between sites, we used the mean hourly temperature averaged across both years (T_avg) and also the maximum temperature recorded in the sampling period, averaged across both years (T_max), because extreme weather events may affect arthropod physiology, behavior, and population and community dynamics (Parmesan, Root, and Willig 2000). We did not compare minimum temperatures because our sampling period did not extend to the coldest months of the growing season.

Solar radiation is often a better indicator of activity level than air temperature for epigeic arthropods (Høye and Forchhammer 2008). We therefore calculated potential annual direct incident radiation (PDIR) using the formula provided by McCune and Keon (2002). PDIR is an index that represents the amount of solar energy that would reach the ground surface if not obscured by clouds or vegetation and is calculated from slope, aspect, and latitude.

Statistical analyses

We performed separate but identical analyses on pollinators and epigeic arthropods, as described in the following sections. Environmental variables were compared among habitats with an analysis of variance of the mean.

Diversity analyses

For epigeic arthropods, the eight pitfall traps and one litter sample from each two-week sampling period were pooled into a single sample for all diversity analyses. Similarly, bee bowls and vane traps from each sampling period were pooled for pollinators. Trap loss was common at some sites (due to damage by wildlife), and many samples included fewer than the full complement of traps. The resulting difference in collecting effort was corrected for as described in the following.

For each habitat type, we explored the completeness of our inventories by comparing the observed species richness (S) to the Chao 1 estimator of S, an asymptotic estimation of the minimum expected species richness (Chao 1984). To measure diversity, we chose the Shannon diversity index (H′), which applies equal weight to all species, and the inverse Simpson’s diversity index (1/D), which weighs common species more heavily. These were calculated with the SpadeR (Chao et al. 2015) and Vegan (Oksanen et al. 2018) packages in R v3.6.2 (R Core Team 2018).

We used generalized linear mixed models to test for differences in diversity and abundance among habitat types with the lme4 R package (Bates et al. 2015). Pitfall traps do not measure true relative abundance of species because they are biased toward more active species (e.g., Topping and Sunderland 1992) and thus more accurately describe “activity–abundance” or “activity–density,” hereafter referred to simply as “abundance.” Models included habitat type as a fixed variable and collection date as well as study site nested within transect as random variables to avoid pseudoreplication. A Poisson distribution was assumed for species richness and abundance data and a lognormal distribution for the noninteger Chao 1 species richness estimator, Shannon diversity index, and inverse Simpson’s diversity index. The regression was weighted by the number of traps per sample to account for the difference in collection effort.

Community analyses

The fifteen sample units (five transects × three habitat types) used in all community analyses represented the pooled abundance of each species within a site over all sampling periods. To reduce the noise in the community data set, we excluded rare species that occurred in less than 10 percent of the sample units (i.e., species occurring at one site only). Species abundances were log-transformed to reduce variation. Outlier species and sites (less than two standard deviations above the mean distance between any two species or sites) were retained if deleting them had no noticeable effect on the analyses. Taxa used in the analyses were species or morphospecies, with the exception of two solitary bee genera (Lasioglossum, Osmia) and one syrphid genus (Eupeodes) that either could not be identified to species or were so rare in our collections that they would have been deleted if species were considered separately. We deemed these taxa ecologically important to consider, with similar enough natural histories within the genus to lump into meaningful taxa for community analyses.

We performed all multivariate community analyses using PC-ORD (v6.22; McCune and Mefford 2011). We looked for patterns in community composition across habitat types using nonmetric multidimensional scaling (NMS; Kruskal 1964) to ordinate fifteen sites in species space and then correlated environmental, climatic, and biological variables of interest (Table 7) with axes in the ordinations. Biological variables included the proportion of social, solitary, or parasitic bee species at a site; functional guilds for syrphid fly larvae; and functional guilds for all epigeic arthropods. NMS is a robust and effective ordination method for nonnormally distributed ecological data with a high proportion of zero values (McCune and Grace 2002). The method uses an iterative search for the best position of N entities (samples or species) on a reduced number of k dimensions (axes) that minimizes the “stress” of the resultant k-dimensional ordination. An ordination with low stress (i.e., less than ten; McCune and Grace 2002) has maximal correspondence between ranked distances in the original n-dimensional space and ranked distances in the reduced k-dimensional space. We used the relative Sørensen distance measure for measuring dissimilarity, which relativizes data points by sample unit total (i.e., absolute abundances for species become proportions of the total pollinator or epigeic arthropod catch at that site), thus helping to account for any difference in sampling efforts between sites. We used habitat type overlays to help interpret patterns of community composition among sites in the ordination plots.

Table 3. The correlation (Pearson’s r) of environmental and biological variables with the two axes displayed in the NMS ordinations of fifteen sites in Denali National Park and Preserve (Figure 4).

Variable	Pollinators		Epigeic arthropods
	Axis 1	Axis 2	Axis 1	Axis 2
Elevation	0.854	−0.005	0.822	0.446
Latitude	−0.238	−0.683	−0.037	−0.126
Longitude	−0.113	−0.878	0.015	0.143
Slope	0.432	0.184	0.504	−0.211
Aspect	0.193	0.180	0.132	−0.240
PDIR	0.583	−0.013	0.654	−0.248
PC1	−0.629	0.144	−0.612	0.149
PC2	0.082	0.252	−0.278	0.146
PC3	0.013	0.200	−0.197	0.062
Mean air temperature (Tavg)	−0.705	0.098	−0.712	−0.598
Max air temperature (Tmax) % Bees total	−0.794 0.864	−0.040 0.627	−0.835	−0.612
% Parasitic bees	0.464	−0.377
% Social bees	0.026	−0.248
% Solitary bees	−0.173	0.363
% Syrphidae total	−0.864	−0.627
% Carnivores (larvae)	0.359	0.437
% Detritivores (larvae)	−0.178	0.068
% Herbivores (larvae)	−0.348	−0.486
% Spiders			−0.309	0.524
% Beetles			−0.014	−0.462
% Grasshoppers			0.593	−0.120
% Detritivores			0.569	0.419
% Herbivores			0.654	−0.107
% Fungivores			0.372	0.103
% Carnivores			−0.818	−0.167

Note. Strong correlations (r > |0.500|) are in bold.

We used blocked multiple-response permutation procedures (MRBP; Biondini, Mielke, and Berry 1988) to test for community composition differences between predefined groups (i.e., forest, shrub, tundra). This nonparametric procedure provides a measure of within-group homogeneity (A = 1 means that all entities within the group are identical) and a p value associated with the null hypothesis of no difference between groups. We used a Euclidean distance measure and medians were aligned to zero for all blocks (transects) to focus the analysis on within-transect differences among habitat types (McCune and Grace 2002).

As a complement to MRBP, we used blocked indicator species analysis (Dufrêne and Legendre 1997) to examine associations of individual species with forest, shrub, and tundra habitat types. This method calculates an indicator value (IV, ranging from 0 to 100) for each species in each habitat type, based on the relative abundance and frequency of the species per habitat type (McCune and Grace 2002). For example, a perfect indicator species for tundra (IV = 100) would be one that is restricted to tundra and occurs in every tundra site. The statistical significance of the IV_max across habitat types for each species was evaluated by a Monte Carlo test using 1,000 randomizations, with a null hypothesis that IV_max is no larger than would be expected by chance alone.

Results

Environmental characteristics

Temperatures and PDIR varied among sites (Table 4). T_avg and T_max were closely correlated. T_avg ranged from 10.4°C at site 2T (tundra) to 13.6°C at site 5S (shrub); T_max ranged from 23.1°C at site 4T (tundra) to 30.3°C at site 5F (forest; Table 4). PDIR ranged from −0.74 at site 4F (forest) to −0.35 at site 1T (tundra; Table 4). T_avg and T_max in tundra were significantly colder than in forest or shrub habitat types, and forest received significantly less PDIR than shrub or tundra habitat types (Figure 3).

Table 4. Environmental characteristics of study sites in Denali National Park and Preserve.

Plot	Slope (°)	Aspect (°)	PDIR	Tavg (°C)	Tmax (°C)	PC1	PC2	PC3
1F	0.3	180	−0.63	13.6 ± 0.1	29.2 ± 1.1	0.07	−0.25	−0.13
1S	15.7	140	−0.42	12.6 ± 0.1	29.2 ± 1.4	−0.26	−0.15	−0.19
1T	20.0	214	−0.35	10.7 ± 0.1	24.7 ± 1.1	−0.28	−0.27	−0.19
2F	13.3	81	−0.68	12.5 ± 0.1	28.1 ± 0.4	0.40	−0.08	0.08
2S	8.3	225	−0.52	12.7 ± 0.1	26.5 ± 1.3	−0.18	0.55	0.29
2T	17.3	116	−0.49	10.4 ± 0.1	23.2 ± 1.6	−0.28	−0.04	0.64
3F	0.0	–	−0.64	12.5 ± 0.1	27.4 ± 0.9	0.37	−0.06	0.03
3S	6.3	214	−0.54	12.1 ± 0.1	26.0 ± 1.8	0.17	0.45	−0.29
3T	26.0	114	−0.46	11.1 ± 0.1	23.0 ± 1.3	−0.24	−0.37	−0.08
4F	5.0	30	−0.74	12.2 ± 0.1	26.6 ± 0.8	0.36	0.04	−0.02
4S	2.7	114	−0.61	12.0 ± 0.1	25.9 ± 1.5	−0.30	0.35	−0.48
4T	6.7	134	−0.54	11.1 ± 0.1	23.1 ± 1.5	0.05	−0.25	0.01
5F	4.7	230	−0.57	13.4 ± 0.1	30.3 ± 0.8	0.33	−0.03	0.02
5S	27.7	99	−0.58	13.6 ± 0.1	29.0 ± 1.0	−0.18	0.01	0.03
5T	10.0	197	−0.47	12.3 ± 0.1	25.3 ± 0.4	−0.03	0.09	0.28

Note. Temperatures are mean daily values ± SE (T_avg) and the average of the highest recorded temperature in 2018 and 2019 (T_max). PDIR was calculated with the formula of McCune and Keon (2002) based on latitude, slope, and aspect. PCs are combinations of fifteen ground cover types from Table 2.

Figure 3. Analysis of variance results for differences in mean and maximum air temperature across forest, shrub, and tundra habitat types in Denali National Park and Preserve. Habitat types with different letters (A, B) indicate significantly different means at α = .05.

The first three principal components (PCs) indicated the overall ground cover composition at each study site and explained 80 percent of the variation among sites (Table 2). PC1, a gradient from spruce forest to open shrubs and tundra, explained 49 percent of the variation. PC2, a gradient of large shrubs to smaller or fewer shrubs, explained 17 percent, and PC3, a gradient from prostrate shrubs to an open moss and lichen understory, explained 13 percent.

Arthropod abundance, richness, and diversity

Our final data set comprised 35,473 adult arthropods, representing 510 species in forty-eight families and five orders (Table 5). We were unable to identify individuals from several families or genera to morphospecies due to taxonomic challenges or physical damage: ptiliid beetles (83 individuals; 100 percent of specimens in this taxon); most aleocharine staphylinid beetles (4,967; 79 percent); Syrphidae (84; 4.6 percent); the syrphid genera Chrysotoxum (4; 100 percent), Eupeodes (113; 92.6 percent), Parasyrphus (5; 2.2 percent), Platycheirus (287; 74.9 percent), Sericomyia (8; 72.7 percent), Sphaerophoria (131; 100 percent), and Sphegina (2; 100 percent); and the bee genera Bombus (101; 1.4 percent) and Lasioglossum (176; 55.7 percent). An additional 125 beetles (0.8 percent) and 59 spiders (0.7 percent) could not be identified beyond order. These specimens were thus included in abundance measures but excluded from diversity estimates.

Table 5. Taxonomic summary of family-level arthropod abundance and species richness across forest, shrub, and tundra habitat types in Denali National Park and Preserve.

	Abundance				Species richness
Taxon	Total	Forest	Shrub	Tundra	Total	Forest	Shrub	Tundra
All taxa	35,473	9,308	12,285	13,880	510	284	333	343
Spiders	8,867	2,519	3,253	3,095	199	108	129	128
Agelenidae	5	2	3	0	1	1	1	0
Amaurobiidae	76	26	19	31	1	1	1	1
Araneidae	11	0	0	11	1	0	0	1
Clubionidae	17	9	5	3	4	4	3	2
Dictynidae	98	17	50	31	7	3	4	2
Gnaphosidae	659	46	396	217	13	4	11	12
Hahniidae	91	49	41	1	2	2	2	1
Linyphiidae	3,793	1,872	968	953	116	75	71	73
Liocranidae	2	0	2	0	1	0	1	0
Lycosidae	3,527	450	1,548	1,529	22	10	15	16
Philodromidae	60	0	14	46	4	0	4	2
Salticidae	32	0	14	18	8	0	6	5
Theridiidae	17	3	1	13	6	3	1	4
Thomisidae	351	26	96	229	11	5	7	9
Titanoecidae	69	0	69	0	2	0	2	0
Unidentified	59	19	27	13	-	-	-	-
Beetles	16,153	5,203	5,984	4,966	210	116	129	145
Brentidae	58	0	0	58	1	0	0	1
Buprestidae	1	1	0	0	1	1	0	0
Byrrhidae	28	13	5	10	7	4	4	5
Byturidae	1	0	0	1	1	0	0	1
Cantharidae	21	1	12	8	1	1	1	1
Carabidae	1,595	804	407	384	32	15	21	22
Cerambycidae	1	0	0	1	1	0	0	1
Chrysomelidae	8	1	1	6	4	1	1	3
Ciidae	4	0	4	0	1	0	1	0
Coccinellidae	22	1	7	14	6	0	3	5
Cryptophagidae	161	74	51	36	12	8	8	5
Curculionidae	297	90	62	145	15	10	7	8
Dytiscidae	1	1	0	0	1	1	0	0
Elateridae	80	22	9	49	9	4	4	8
Eucinetidae	7	0	7	0	1	0	1	0
Hydrophilidae	3	0	0	3	1	0	0	1
Lampyridae	3	1	0	2	-	-	-	-
Latridiidae	332	24	118	190	12	5	9	10
Leiodidae	224	35	116	73	11	5	7	10
Mordellidae	2	0	1	1	2	0	1	1
Nitidulidae	20	3	3	14	3	3	1	1
Ptiliidae	83	43	34	6	-	-	-	-
Ptinidae	12	1	2	9	1	1	1	1
Scarabaeidae	486	7	1	478	2	1	1	2
Silphidae	124	9	49	66	2	1	2	1
Staphylinidae	12,454	4,020	5,042	3,392	83	55	56	58
Unidentified	125	52	53	20	-	-	-	-
Syrphid flies	1,836	674	503	659	67	43	48	44
Syrphidae	1,836	674	503	659	67	43	48	44
Bees	7,734	912	2,491	4,331	30	17	24	22
Andrenidae	7	2	4	1	1	1	1	1
Apidae	7,375	863	2,193	4,319	18	12	15	17
Colletidae	25	13	12	0	1	1	1	0
Halictidae	316	33	275	8	4	2	4	1
Megachilidae	11	1	7	3	6	1	4	3
Grasshoppers	883	0	54	829	4	0	3	4
Acrididae	883	0	54	829	4	0	3	4

Bumble bees (Bombus) accounted for 77 percent of total pollinator abundance and 95 percent of all bees collected. They included the seven most abundant pollinator species, and two species, Bombus lapponicus sylvicola and B. mixtus, respectively accounted for 25 and 16 percent of all pollinators. Syrphid flies were represented by more than twice as many species as bees, even though they made up less than 20 percent of the total pollinator catch (Table 5). Two species, Melanostoma mellinum and Meliscaeva cinctella, together made up 23 percent of the total syrphid catch. The epigeic arthropod community was more even, because no species accounted for more than 5 percent of individuals. The most abundant species included wolf spiders (Pardosa lapponica, 5 percent; Pardosa uintana, 3 percent), beetles (Eucnecosum brunnescens, 5 percent; Pterostichus brevicornis, 3 percent), and a grasshopper (Aeropedellus arcticus, 3 percent). Roughly one-quarter of all arthropod species (132) were represented by a single individual. Abundance did not differ significantly among habitat types for epigeic arthropods or syrphid flies, but significantly more bees were collected in shrub and tundra than in forest (Table 6).

Table 6. Mean abundance, observed species richness (S), estimated species richness (Chao 1), Shannon diversity index (H′), and inverse Simpson’s index of arthropods per sample (one collection at one site) in forest, shrub, and tundra habitat types in Denali National Park and Preserve.

	Forest	Shrub	Tundra	χ^2	df	p
Spiders
Abundance	48.4 ± 5.4	62.6 ± 6.2	59.5 ± 7.3	1.86	2	.3945
Observed richness	14.6 ± 0.9	18.1 ± 1.1	15.1 ± 1.0	4.64	2	.0983
Estimated richness	23.9 ± 1.7	28.8 ± 1.6	24.6 ± 1.9	3.09	2	.2135
Observed/estimated richness	61.1%	62.8%	61.4%
Shannon diversity	1.1 ± 0.1	2.3 ± 0.1	2.1 ± 0.1	4.30	2	.1165
Inverse Simpson’s diversity	7.0 ± 0.5	8.3 ± 0.5	5.5 ± 0.3	5.99	2	.0501
Beetles
Abundance	100.1 ± 21.5	115.1 ± 22.6	95.5 ± 10.8	0.40	2	.8176
Observed richness	13.9 ± 1.2	15.5 ± 1.2	17.2 ± 1.1	1.41	2	.4933
Estimated richness	23.4 ± 2.2	25.0 ± 2.1	30.2 ± 2.7	0.89	2	.6400
Observed/estimated richness	59.4%	62.0%	57.0%
Shannon diversity	1.9 ± 0.1	2.1 ± 0.1	2.2 ± 0.1	2.16	2	.3397
Inverse Simpson’s diversity	5.6 ± 0.4	7.2 ± 0.5	7.3 ± 0.6	3.26	2	.1960
Syrphid flies
Abundance	14.0 ± 2.7	10.5 ± 1.4	13.7 ± 2.8	0.74	2	.6913
Observed richness	4.1 ± 0.5	3.6 ± 0.4	2.6 ± 1.0	5.48	2	.0645
Estimated richness	5.7 ± 0.8	6.4 ± 1.0	4.8 ±1.0	0.39	2	.8242
Observed/estimated richness	71.9%	56.3%	54.2%
Shannon diversity	0.9 ± 0.1	0.9 ± 0.1	0.6 ± 0.1	076	2	.6850
Inverse Simpson’s diversity	2.9 ± 0.3	3.1 ± 0.2	2.4 ± 0.3	5.41	2	.0670
Bees
Abundance	19.0 ± 3.1^A	51.9 ± 7.5^B	90.2 ± 11.4^B	15.26	2	.0005*
Observed richness	3.5 ± 0.4^A	6.6 ± 0.5^B	8.2 ± 0.3^B	47.52	2	<.0001*
Estimated richness	4.1 ± 0.5	7.8 ± 0.6	9.4 ± 0.5	2.76	2	.2512
Observed/estimated richness	85.4%	84.6%	87.2%
Shannon diversity	0.8 ± 0.1	1.4 ± 0.1	1.5 ± 0.1	2.96	2	.2597
Inverse Simpson’s diversity	2.3 ± 0.2^A	3.7 ± 0.2^B	3.6 ± 0.2^B	7.96	2	.0187*
Grasshoppers
Abundance	0	1.0 ± 0.5	15.9 ± 5.1	3.92	1	.1412
Observed richness	0	0.2 ± 0.1	0.6 ± 0.1	1.77	1	.4119
Estimated richness	0^AB	0.3 ± 0.1^B	0.6 ± 0.1^BC	6.77	1	.0340*
Observed/estimated richness	—	66.7%	100%
Shannon diversity	0	0.0 ± 0.0	0.1 ± 0.2	1.47	1	.4792
Inverse Simpson’s diversity	0	1.2 ± 0.1	1.2 ± 0.1	0.01	1	.9273

Notes. p Values based on generalized linear mixed models.

*Indicates significant differences at α = 0.05. Capital letters indicate significantly different means.

Table 7. Pollinator and epigeic arthropod taxa from Denali National Park and Preserve with strong correlations (r > |0.500|) with axes from the NMS ordination (Figure 4) and strong indicator species (IV> 33, p < .05) of forest, shrub, and tundra habitat types from blocked indicator species analysis.

Taxon	NMS		Indicator species analysis
	Axis 1	Axis 2	Indicator group	Indicator value	p Value
Syrphid flies
Chalcosyrphus vecors		−0.565
Cheilosia latrans	−0.615	−0.575
Chrysotoxum sp. 1	−0.592
Dasysyrphus intrudens		−0.613
Dasysyrphus nigricornis	−0.600	−0.654
Dasysyrphus venustus	−0.526	−0.635
Eupeodes spp.		0.612	Tundra	60.1	.038
Lapposyrphus lapponicus		−0.574
Megasyrphus laxus	−0.593	−0.744
Melanostoma mellinum	−0.775
Meliscaeva cinctella	−0.855	−0.571	Forest	69.6	.011
Parasyrphus macularis	−0.747		Forest	75.4	.035
Parasyrphus relictus	−0.636
Parasyrphus vockerothi	−0.627
Platycheirus sp. 1	−0.502
Sericomyia arctica/jakutica	−0.633
Syrphus torvus		−0.596
Xylota subfasciata		−0.505
Bees
Bombus cryptarum	0.583
Bombus flavidus	0.712
Bombus flavifrons	0.786
Bombus frigidus	0.833
Bombus insularis	0.685
Bombus jonellus		0.808
Bombus kirbiellus	0.913		Tundra	60.2	.015
Bombus kluanensis	0.542
Bombus lapponicus sylvicola	0.750	0.600	Tundra	45.0	.015
Bombus melanopygus		−0.782
Bombus natvigi	0.642		Tundra	80.0	.038
Bombus neoboreus	0.759		Tundra	69.7	.015
Bombus polaris	0.774		Tundra	79.8	.015
Osmia spp.		0.595
Spiders
Acanthoneta aggressa	−0.555
Agelenopsis utahana		−0.579
Agyneta allosubtilis			Tundra	77.7	.036
Agyneta jacksoni	0.673		Tundra	80.0	.038
Arctella lapponica	0.781	0.539	Tundra	69.3	.035
Bathyphantes eumenis	−0.522		Forest	56.8	.037
Bathyphantes pallidus		−0.930
Ceraticelus fissiceps	−0.578		Forest	73.2	.039
Chalcoscirtus alpicola	0.709
Clubiona bryantae		−0.541
Clubiona kulczynskii	−0.601
Cnephalocotes obscurus			Shrub	84.0	.030
Diplocentria bidentata	−0.595
Diplocentria rectangulata	−0.522
Diplocephalus subrostratus	−0.516
Enoplognatha intrepida	0.561
Evarcha proszynskii		−0.568
Gnaphosa borea	0.695
Gnaphosa microps			Shrub	50.0	.017
Gnaphosa orites	0.774		Tundra	80.0	.035
Gonatium crassipalpum	0.631
Hahnia glacialis	−0.508
Hilaira herniosa	−0.741		Forest	57.5	.024
Horcotes strandi	0.750		Tundra	69.9	.041
Improphantes complicatus	−0.766
Islandiana cristata		0.601
Ivielum sibiricum			Shrub	58.3	.021
Lepthyphantes alpinus	−0.702		Forest	47.7	.045
Linyphiidae DNP 01	0.719
Micaria aenea		−0.649
Micaria alpina	0.854		Tundra	71.6	.039
Micaria constricta		−0.526
Microneta viaria		−0.519
Ozyptila arctica	0.635
Ozyptila sincera		−0.524	Shrub	69.0	.014
Pardosa lapponica	0.779		Tundra	62.1	.027
Pardosa moesta			Shrub	58.9	.017
Pardosa podhorskii	0.864		Tundra	69.9	.035
Pardosa uintana	−0.820
Sciastes truncatus	−0.627
Scotinotylus sacer	−0.575
Sisicottus montanus	−0.743		Forest	77.7	.013
Sisicottus nesides		−0.62
Sisicottus orites	−0.541
Styloctetor compar		−0.502
Thanatus arcticus	0.794		Tundra	70.3	.050
Tiso aestivus	0.646
Titanoeca nivalis		−0.653
Walckenaeria atrotibialis	−0.685	−0.678
Walckenaeria castanea		−0.578
Walckenaeria directa		−0.714
Xysticus emertoni	0.744		Tundra	64.3	.013
Xysticus keyserlingi		−0.571
Xysticus luctuosus		−0.618
Zelotes sula		−0.714
Zornella armata	−0.827	−0.557
Beetles
Acidota quadrata	0.515
Agoliinus piceatus	0.712
Aleochara castaneipennis		−0.754
Amara alpina	0.780		Tundra	100.0	.013
Amara brunnea			Shrub	89.8	.017
Amara hyperborea		−0.727
Ascoliocerus sanborni	0.567		Tundra	80.0	.038
Atheta altaica		−0.548
Atheta munsteri	0.655
Atheta prudhoensis	0.654
Atheta strigosula		−0.639
Atomaria kamtschatica			Shrub	80	.045
Bolitobius cingulatus		−0.641
Byrrhus kirbyi		−0.515
Caenocara scymnoides	0.682
Carabus chamissonis	0.678		Tundra	83.1	.024
Cartodere constricta		−0.505
Catops alpinus			Shrub	64.1	.028
Ceutorhynchus DNP 1	0.655
Chaetocnema protensa	0.547
Corticaria DNP 1	0.702		Tundra	63.5	.009
Corticaria DNP 2	0.628
Corticaria DNP 3		−0.651
Cryptophagus difficilis		−0.529
Denticollis varians		−0.547
Dichelotarsus sp		−0.515
Hypera diversipunctata	0.521		Tundra	71.7	.038
Ischnosoma fimbriatum		−0.656
Leiodes rhaetica	0.526
Liogluta granulosa		−0.526
Loborhynchapion cyanitinctum	0.768		Tundra	42.8	.038
Lordithon fungicola		−0.718
Lypoglossa angularis		−0.61
Megarthrus angulicollis	−0.626
Micropeplus laticollis	−0.505	−0.613
Neothetalia pallitarsis		−0.53
Olophrum latum		−0.847
Omalium DNP 3	0.552		Tundra	70.0	.039
Oxypoda sylvia	−0.585	−0.636
Proteinus parvulus	−0.646	−0.732
Pterostichus barryorum		−0.567
Pterostichus pinguedineus	−0.511
Pterostichus riparius		−0.578
Quedius brunnipennis	−0.627
Quedius fellmani		0.603
Quedius labradorensis		−0.55
Scaphium castanipes		−0.637
Scymninae DNP 1	0.621
Sericoda quadripunctata		−0.571
Sericus incongruus		−0.515
Silusa californica		−0.62
Stenichnus DNP 1		−0.631
Stenus austini	−0.687		Forest	72.1	.035
Sylvanelater mendax		−0.598
Tachinus basalis	−0.506
Tachinus elongatus			Tundra	51.9	.013
Tachinus frigidus		−0.753
Tachinus instabilis	0.614
Tachinus quebecensis		−0.856
Tachyporus nimbicola	0.513
Tachyporus rulomus	0.615
Grasshoppers
Aeropedellus arcticus	0.773		Tundra	71.4	.027
Bohemanella frigida		−0.501

Across habitat types, tundra had the most unique species with ninety-five, followed by shrubs with sixty-nine and forest with fifty-one. Observed species richness and diversity (inverse Simpson’s) were significantly lower in forest than in shrub or tundra for bees, and although not significant, there was a strong trend of decreasing species richness for syrphid flies from forest to tundra (Table 6). Among epigeic arthropods, the estimated species richness of grasshoppers was significantly lower in forest than in tundra. Beetle diversity showed a nonsignificant trend of increasing from forest to tundra, whereas spider diversity peaked in shrubs (Table 6).

Community analyses

Pollinators

We started with a matrix of ninety-three pollinator species (sixty-nine syrphids, twenty-four bees) × fifteen sites for all multivariate community analyses and then deleted twenty rare species (sixteen syrphids, four bees), resulting in a matrix of seventy-three species (these included the lumped genera Eupeodes, Andrena, Osmia, and Lasioglossum) × fifteen sites. The NMS two-dimensional solution represented 92.5 percent of the variation in the original data and had a stress value of 8.4 (Monte Carlo test, p = .004), achieved after forty-four iterations. We rotated the ordination to maximize differences between treatments (habitat types) along axis 1, which resulted in 49.7 percent of the original variation being loaded onto axis 1 and 42.8 percent onto axis 2 (Figure 4a).

Figure 4. Nonmetric multidimensional scaling ordination of fifteen sites in species space for (a) pollinators and (b) epigeic arthropods in Denali National Park and Preserve. Habitat type overlays are color-coded.

The three habitat types were separated out along axis 1, with some overlap between sites in tundra and shrub habitat types (Figure 4a). Strong positive correlations (r > 0.5) with axis 1 were seen for elevation and PDIR, whereas mean and maximum air temperature and the ground cover principal component PC1 were strongly negatively correlated (r < 0.5) with axis 1 (Table 3). Both longitude and latitude were strongly negatively correlated with axis 2 (Table 3).

The proportion of bees at a site increased along axis 1 (i.e., higher proportions in shrub and tundra sites), with a corresponding decrease in the proportion of syrphid flies (Table 3). Among individual taxa, twelve syrphids had a strong negative correlation with axis 1 (i.e., higher proportions in forested sites), ten syrphids had a strong negative correlation with axis 2, and syrphids in the genus Eupeodes had a strong positive correlation with axis 2 (Table 7). For bees, eleven species had strong positive correlations with axis 1, three taxa were strongly positively correlated with axis 2, and Bombus melanopygus was strongly negatively correlated with axis 2 (Table 7).

The separation of habitat types shown in the NMS plot was supported by MRBP, which revealed a significant difference in species composition between habitat types (within-group homogeneity, A = 0.203, p < .001). Pairwise comparisons showed significant differences between all habitat types: forest vs. shrub (A = 0.168, p = .016), forest vs. tundra (A = 0.314, p = .015), and shrub vs. tundra (A = 0.195, p = .016).

Indicator species analysis revealed six strong indicator taxa for tundra, including Eupeodes syrphid flies and five bumble bee species (Table 7). In addition, there were two forest indicators among the syrphid flies, Meliscaeva cinctella and Parasyrphus macularis (Table 7).

Epigeic arthropods

We started with a matrix of 412 species (199 spiders, 209 beetles, 4 grasshoppers) × fifteen sites for all multivariate community analyses and then deleted 147 rare species (76 spiders, 70 beetles, 1 grasshopper), resulting in a matrix of 265 species × fifteen sites. The NMS two-dimensional solution represented 90.0 percent of the variation in the original data and had a stress value of 9.2 (Monte Carlo test, p = .004), achieved after seventy-seven iterations. We rotated the ordination to maximize differences between treatments (habitat types) along axis 1, which resulted in 78 percent of the original variation being loaded onto axis 1 and 12.2 percent onto axis 2.

The three habitat types were clearly separated out along axis 1 (Figure 4b). Elevation, slope, and PDIR were strongly correlated with axis 1, whereas PC1 and air temperature were strongly negatively correlated with axis 1 (Table 3). Air, soil, and ground temperature were all strongly negatively correlated with axis 2 (Table 3).

The proportions of grasshoppers, detritivores, and herbivores at a site were all strongly positively correlated with axis 1 (i.e., highest proportions in tundra sites), whereas carnivores were strongly negatively correlated with axis 1 (i.e., highest proportions in forest sites; Table 3). The proportion of spiders was strongly positively correlated with axis 2.

Among individual species, thirty-seven were strongly positively correlated with axis 1 (highest proportions in tundra sites: twenty beetles, sixteen spiders, one grasshopper), whereas twenty-six species were strongly negatively correlated with axis 1 (highest proportions in forest sites: eighteen spiders, eight beetles; Table 7). One beetle and two spiders were strongly positively correlated with axis 2, and eighteen spiders were strongly negatively correlated with axis 2 (Table 7).

As we found for pollinators, the strong difference in community composition between the three habitat types shown in the NMS plot was supported by MRBP, which revealed a significant difference in species composition between habitat types (within-group homogeneity, A = 0.141, p < .001). Pairwise comparisons showed significant differences between all habitat types: forest vs. shrub (A = 0.122, p = .015), forest vs. tundra (A = 0.219, p = .015), and shrub vs. tundra (A = 0.154, p = .016).

Indicator species analysis revealed nineteen strong indicators for tundra (ten spiders, eight beetles, one grasshopper), eight strong shrub indicators (five spiders, three beetles), and six strong forest indicators (five spiders, one beetle; Table 7). Amara alpina, a carabid beetle, was the only species with a perfect indicator value (IV = 100), indicating that it was collected in all five tundra sites and in no other habitat types.

Discussion

Arthropod distributions across habitat types and elevation

High latitudes and high elevations are experiencing intense ecological changes in response to anthropogenic climate change (Parmesan 2006; Myers-Smith et al. 2011; Mekonnen et al. 2021). Upslope shifts in vegetation zones have the potential to change arthropod community structure and function. We compared patterns of diversity and community composition for both epigeic and aerial arthropods across an elevation gradient in DNPP and examined their relationships with vegetation cover variables.

Our study revealed that species composition of both pollinators and epigeic arthropods differed significantly between forest, tall shrub, and alpine tundra habitat types along an elevation gradient in DNPP. Differences were most pronounced for epigeic arthropods, many of which have low dispersal power. Epigeic communities also varied in functional guild structure across habitat types and elevation, with carnivores more dominant in low-elevation forests and herbivores and detritivores making up a higher proportion of tundra communities. Across each of the five taxonomic orders surveyed, the distributions of a large proportion of species (23–55 percent) were strongly correlated with the elevation gradient.

Among pollinators, bees dominated most shrub and tundra sites in species richness and abundance, and bumble bees (Bombus) made up 95 percent of our bee collections. Though a few generalist species were fairly common across all habitat types (e.g., Bombus cryptarum, B. mixtus), a suite of bumble bee species were strong indicators of tundra. These included four species in the circumpolar subgenus Alpinobombus (B. kirbiellus, B. natvigi, B. neoboreus, B. polaris) and also B. lapponicus sylvicola, the most abundant bee in our study. Bees were the only arthropod group to show significantly higher abundance, species richness, and diversity in shrub and tundra versus forest.

The overlap in species composition between some shrub and tundra sites suggests that some bumble bee species may be using both habitat types for foraging. Given the compressed spatial distances between low-, mid-, and high-elevation vegetation communities in these northern mountains relative to the foraging distances that bees are able travel (Greenleaf et al. 2007; Zurbuchen et al. 2010), it may be, for example, that overwintering queen bumble bees and the first broods of workers forage on early blooming tall willow flowers in the spring (Vogt and Heinrich 1994) whereas later in the season workers are primarily foraging on forbs and dwarf shrubs in the tundra, where vascular plant species richness is highest (Roland and Schmidt 2015). Our sampling methods for pollinators did not allow us to track seasonal changes in species composition.

Solitary bees made up less than 5 percent of the total bee catch and were most abundant and diverse in shrub sites. Most of these bees are ground nesters (Andrena, Anthophora, Lasioglossum, Osmia), a few nest in dead wood or hollow stems (Hylaeus, Megachile, Osmia), and we collected one cleptoparasite (Coelioxys funeraria) of leafcutter (Megachile) bees. Solitary bees were patchy in distribution and thus not very informative for our community analyses: one shrub site accounted for 66 percent of the total solitary bee abundance, including 66 percent of all of the solitary taxa collected. Soil-nesting sweat bees in the genus Lasioglossum were the most abundant and widespread solitary group, with most collected in shrub sites and only a single specimen from the tundra. Masked bees (Hylaeus annulatus), which nest in hollow stems, occurred exclusively in forest and shrub sites. Previous bee surveys in interior and arctic Alaska (Armbruster and Guinn 1989; Rykken 2015, 2017) have revealed that many solitary and cleptoparasitic bees are very patchy in their distributions, favoring disturbed and open habitat types not represented by our transects (e.g., sandy river bluffs, floodplains, roadsides); thus, it was not surprising that solitary bees were not strong indicators of any of the habitat types examined in this study.

Syrphid flies made up a higher proportion of the pollinator diversity in forest sites than in shrub or tundra, and two species, Meliscaeva cinctella and Parasyrphus macularis, were identified as strong indicators of forest. As a general rule, dipteran diversity and abundance are known to increase with latitude and elevation (Kearns 1992; Loboda et al. 2017; McCabe and Cobb 2021). Focusing only on syrphid flies in this study limited our view of overall dipteran diversity and abundance but allowed us to use species-level identifications in our analyses, which can be a challenge for many dipteran families. Field observations and by-catch in the samples suggested that muscoid flies far outnumbered syrphid flies across all habitat types, but relative proportions of these and other flies were not quantified.

Grasshopper abundance showed a strong elevational trend. They were abundant in tundra, rare in shrubs, and absent from all forests, but the assumptions for statistical tests were not met because they were relatively uncommon and their abundance varied greatly. Only two species, Aeropedellus arcticus and Bohemanella frigida, were represented by more than three specimens, and B. frigida was restricted to transect 5. Aeropedellus arcticus was widespread and among the most abundant arthropods. Both species are geographically restricted to Beringia (Alaska and northwestern Canada; Catling 2008), and Beringian distributions are often associated with relictual tundra areas (Roland and Schmidt 2015; M. E. Edwards, Lloyd, and Armbruster 2018).

We did not observe significant differences in the abundance or diversity of beetles or spiders across the gradient. An investigation of three mountain ranges in Yukon, Canada, also found no consistent effect of elevation on spider richness (Bowden and Buddle 2010b). Previous work along a gradient of shrub dominance in Alaska found that spiders were more abundant in open tundra, but beetle abundance varied by transect, and the overall explanatory power of vegetation height was low (Rich, Gough, and Boelman 2013). Pitfall traps favor ground-active arthropods, such as wolf spiders (Lycosidae). Web-building spiders and xylophagous beetles are abundant in boreal forests but poorly represented in pitfall catches because they spend little time on the ground. Overall beetle and spider abundance and diversity may therefore be higher at forest sites than our pitfall catches suggest.

Spiders showed strong preferences for habitat type, with twenty species found to be significant indicators of single habitat types. The sheet web spiders (Linyphiidae) were the only spiders indicative of forests, where they were the most abundant family by far. Active hunters, like wolf spiders, were generally more associated with shrubs and tundra, perhaps because open ground facilitates prey-finding.

The ordinations for both epigeic arthropods and pollinators revealed strong block effects (i.e., differences in community composition between transects) along axis 2. For pollinators, there was a strong longitudinal effect, with most sites positioned along an east–west gradient. Fifteen pollinator taxa were strongly correlated with axis 2, indicating that for these species, proportions were higher in sites/transects at either the eastern or western end of our study area. None of the environmental variables we measured were strongly correlated with this longitudinal gradient; thus, it is unclear what factors are driving this pattern for pollinators. For epigeic arthropods, sites in transects 1 and 5 were at one end of axis 2, whereas sites in transects 2, 3, and 4 were at the other end. Both mean and maximum air temperature were strongly associated with the gradient along axis 2, suggesting that differences in air temperature across transects are also influencing community composition. Despite these strong secondary block effects, differences in community composition across habitat types and along an elevation gradient emerged clearly for both pollinators and epigeic arthropods.

Environmental drivers

Mean air temperature, maximum air temperature, and PDIR were strongly correlated with pollinator and epigeic arthropod distributions across habitat types and elevation. Viterbi et al. (2013) also found that air temperature was a strong driver for species composition of beetles, butterflies, and spiders along an elevation gradient in the Italian Alps. Temperature may have both direct and indirect effects on structuring arthropod communities. Critical high and low thermal limits are known drivers of species distributions across latitude, with high-latitude species having predictably wider thermal tolerances due to larger variation in seasonal and diurnal temperatures (Sunday et al. 2019). Generally, upper thermal limits are less variable, whereas lower limits decline with increasing latitude (Addo-Bediako, Chown, and Gaston 2000). Critical thermal maxima and minima have been less studied along elevation gradients, but Oyen, Giri, and Dillon (2016), in a relevant study on bumble bee thermal limits in Wyoming, found that high-elevation species (including a species common in our data, Bombus lapponicus sylvicola) were more cold tolerant and less heat tolerant than lower elevation species. However, Slatyer and Schoville (2016) found that thermal limits did not explain turnover of ground beetles along elevation gradients in the Pacific Northwest and surmised that biotic interactions or chronic thermal stress may be setting elevational limits for some species. Behavior and microclimate (e.g., seeking shelter in burrows) may also compensate for elevation-driven declines in temperature (Dillon et al. 2006; Sunday et al. 2014). Thus, though temperature was strongly correlated with insect distributions across our sites, calculating the actual thermal tolerances of species along elevation gradients may be far more complex than measuring their critical maxima and minima.

Even as temperature maxima and minima may set critical thermal limits for insect development and survival, a narrower range of diurnal temperatures is also important for determining activity niches (Guo et al. 2020), including arthropod flight and foraging activities (Taylor 1963; Kenna, Pawar, and Gill 2021). Bumble bees are heterothermic, able to warm their flight muscles on cool days, but ambient temperatures constrain flight performance and optimal thermal niches vary for bees of different sizes and species (Corbett et al. 1993; Kenna, Pawar, and Gill 2021). Other alpine insects, like beetles and syrphid flies, must rely on external sources to warm up their flight muscles but likely also have thermal optima for flight (Taylor 1963; Peng, Fletcher, and Sutton 1992).

Air temperature had a strong influence on within-habitat type variation of epigeic arthropod communities. Within all three habitat types, colder sites had a higher proportional abundance of spiders and fewer beetles than warm sites. Our findings may reflect genuine differences in abundance or could be a sampling artifact caused by lower activity of beetles at cold sites, which would lead to fewer captured in pitfall traps.

Temperature may also affect arthropod communities indirectly through its role in structuring plant communities, and ground cover (PC1) was another strong driver for arthropod distributions in this study. Roland et al. (2019) suggested that microclimate controls the distribution, abundance, and richness of upright woody plants (i.e., forest, tall shrubs) in interior Alaska, but other environmental factors such as ground disturbance and soil pH are stronger drivers of vegetation such as dwarf shrubs, graminoids, and forbs associated with tundra. In DNPP, vascular plant species richness peaks in high-elevation tundra (Roland and Schmidt 2015), providing plentiful nectar and pollen for both syrphid flies and bees. Dense tree and shrub canopy limits the abundance and diversity of forbs available for pollinator foraging but may provide essential resources for some pollinator species such as syrphid flies whose larvae feed on dead wood or sap or mason and leafcutter bees that use preexisting holes made by wood-boring beetles and other insects for nesting. Many syrphid larvae are predaceous on aphids and other soft-bodied plant-feeding insects, thus associated with the host plants of their prey. Web-building spiders use vegetation as attachment points for their webs, and vegetation structure was the primary driver of spider assemblages across a latitudinal gradient in Canada (Bowden and Buddle 2010a).

Woody plants exert additional influence on epigeic arthropod communities through shading (Seibold et al. 2016), generating leaf litter and debris (Uetz 1979; Seibold et al. 2016), accumulating snow, and delaying spring snowmelt (Rich, Gough, and Boelman 2013). Date of snowmelt was negatively correlated with overall arthropod abundance (Asmus et al. 2018) and spider abundance (Rich, Gough, and Boelman 2013) in arctic Alaska. A study in Greenland correlated late snowmelt with smaller body size in spiders (Høye et al. 2009; Bowden et al. 2015), which may lead to reduced fecundity (Ameline et al. 2018; Fox and Czesak 2000) and, in turn, lower abundance. However, snow fence experiments in subarctic Canada found that increased snow depth led to greater body mass and ovary size in female P. lapponica, the dominant spider in our system, likely due to higher subnival temperatures for overwintering spiders (Leguat and Weis 2013). A similar link between growing season length and body mass was proposed for the tundra-dwelling carabid beetle Amara alpina (Becker et al. 2020).

Assessing effects of climate change on arthropod distributions

The years 2016 and 2018 were the warmest on record for central Alaska, and in DNPP, mean annual air temperature increased by more than 2°C between 2006 and 2019 (Sousanes and Hill 2020; Swanson, Sousanes, and Hill 2021). During the last century, as a consequence of climate change, the forest–tundra ecotone in the Alaska Range has advanced 150 m upslope (Lloyd and Fastie 2003) and it is now estimated that shrubs are displacing tundra at a rate of 5 percent per decade in southcentral Alaska (Dial et al. 2007). In DNPP, alpine tundra and its rich flora are spatially restricted compared to lowland forest, because the conditions promoting high plant diversity (cool summer temperatures, ground disturbance, weakly acidic soils) occur in a narrow band of elevation (Roland et al. 2019). The upper limit of the tundra is typically bounded by steep, rocky, unstable slopes with little soil development and thus upslope colonization by tundra plants will likely progress at a much slower rate than the upslope migration of trees and shrubs. Thus, of the three habitat types we sampled in this study, tundra stands out as the most vulnerable to effects from climate change.

In a study area that overlaps ours in DNPP, Mizel et al. (2016) reported that passerine birds changed their elevation distribution over a period of eighteen years, with pervasive upward shifts for species associated with the shrub–tundra ecotone and somewhat weaker responses from forest-associated species. The authors surmised that because shrub colonization into the tundra is happening at a faster rate than tree recruitment upslope, birds associated with these different ecotones are responding at comparable rates. As climate trends continue, we predict a corresponding change in arthropod communities as some species track their preferred climate niche or habitat boundaries. Some forest-associated species may expand their ranges, and some shrub-dwelling arthropods may also push upslope. Loss of connectivity between alpine areas due to encroaching shrubs and trees may isolate some tundra populations (Barredo, Mauri, and Caudullo 2020; Inouye 2020). Different response rates or patterns across species could alter interactions and change community structure and function (Petchey et al. 1999). This includes mismatches in phenology between plants and their pollinators if they respond to earlier snowmelt and warming at different rates (Kudo 2014; Pyke et al. 2016). Species loss in the alpine could have amplified consequences because tundra communities are thought to have low functional redundancy (Petchey et al. 1999).

In our study, tundra had more unique species and strong indicator species than any other habitat type. Among bees, four of the strong tundra indicators were bumble bees in the subgenus Alpinobombus. This group of arctic bees is circumpolar in distribution, with five species known from North America and four from Eurasia (Williams et al. 2019). Rasmont et al. (2015) modeled future distributions of the four European Alpinobombus species under several climate change scenarios of increasing severity and concluded that all were at very high climate change risk, with projections of extirpation or extinction for three of the species by 2100. Lee, Williams, and Pearson (2019) assessed climate vulnerability for North American Alpinobombus species and concluded that they were at a higher risk of decline than more temperate alpine bumble bee species, illuminating the fact that in DNPP, some high-elevation tundra species are at increased risk because their range is also restricted to northern latitudes. Notably, a fifth Alpinobombus species in DNPP, Bombus kluanensis, was also restricted to tundra but too rarely captured to meet the criteria for an indicator species in this study. The species was recently described (Williams and Cannings 2016) and its known range—among the smallest of any North American Bombus species—DNPP in Alaska and the Kluane region of the Yukon, Canada—suggests a Beringian distribution.

In contrast, the other strong alpine tundra indicator among bees in our study was Bombus lapponicus sylvicola. This species comprised 30 percent of all bees collected and was found at every site. It occurred rarely in forests, moderately in shrubs, but most abundantly in alpine tundra. Bombus lapponicus sylvicola is found across Alaska and Canada and ranges south into the Sierra and Rocky Mountains, where it occurs at higher elevations. Its high abundance relative to any of the Alpinobombus species and its distribution across the three habitat types indicate that it may be far more adaptable to shifting vegetation boundaries caused by climate change.

Obligate tundra species are most likely to experience future range contractions (Barredo, Mauri, and Caudullo 2020), although this is dependent on mobility and dispersal ability. The carabid beetle, Amara alpina, the only perfect tundra indicator in our study, is a Holarctic ground beetle typically found in tundra. It is believed not to fly, even though adults have fully formed wings (Lindroth 1969), and its mobility is further constrained by the nonvagile egg stage, unlike some spiders that carry their egg sacs. Amara alpina may therefore be less capable of tracking suitable environmental conditions or maintaining gene flow between habitat patches (Beckers et al. 2020). This hypothesis is supported by genetic evidence that revealed that the rapid contraction of tundra at the end of the Pleistocene led to the extinction of some populations and isolation of others on alpine peaks in the Appalachian and Rocky Mountains (Reiss, Ashworth, and Schwert 1999). Another tundra indicator beetle, Tachinus elongatus, shares a similar North American distribution to A. alpina (Campbell 1973), probably due to the same historical processes. These beetles may be indicators for a suite of low-vagility tundra specialists whose ranges are likely to further fragment and shrink due to anthropogenic climate change. In contrast, wolf spiders (e.g., Pardosa spp.) are highly mobile in all life stages and capable of dispersing long distances via ballooning (Yip and Rayor 2014). They may be more resilient to environmental change because they can quickly move into favorable microhabitats and may even increase in abundance or expand their distribution (Becker et al. 2020).

A call for long-term monitoring

Historical baseline data that can be used for comparisons with present-day arthropod distributions along elevation gradients are rare (Franzén and Öckinger 2012; Cerrato et al. 2019; Halsch et al. 2021), pointing to the urgent need for long-term monitoring of arthropods in alpine systems. Arthropods are an essential and diverse component of all but the highest montane environments, and disruptions to their communities will have consequences for fundamental ecosystem processes such as pollination, decomposition, predator–prey dynamics, and herbivory. In a study similar to ours in the Italian Alps, Viterbi et al. (2013) found that alpine zones above tree line harbored more endemic and vulnerable pollinator and epigeic arthropod species than any of the lower zones they sampled, thus highlighting the need for a conservation and monitoring program in this system. There is a growing consensus among entomologists working in montane systems worldwide that long-term monitoring is urgently needed (Franzén and Öckinger 2012; Meyer et al. 2015; Inouye 2020). Ideally, a global network of monitoring sites such as those used by programs currently focused on alpine plants (Global Observation Research Initiative in Alpine environments (GLORIA)); Grabherr, Gottfried, and Paull 2000) or arctic tundra arthropods (Network for Arthropods of the Tundra (NeAT)); Høye and Culler 2018) would capture trends for different arthropod taxa across latitudinal as well as elevation gradients.

For individual monitoring programs to be effective and informative, it is critical to include arthropod taxa that are diverse and abundant across the elevation gradient. Our study in subarctic Alaska included several taxa that are often overlooked due to taxonomic challenges but appeared to be strong indicators of individual habitat types. For example, pollinator studies along elevation gradients have tended to focus on bees and butterflies (Franzén and Öckinger 2012; Pyke et al. 2016; Cerrato et al. 2019) in part because they are conspicuous and relatively well known taxonomically. However, syrphid flies, bee flies (Bombyliidae), dance flies (Empididae), muscoid flies (including Muscidae, Anthomyiidae), and others become an increasingly abundant and important component of the pollinator fauna at higher latitudes and elevations (Kearns 1992; Elberling and Olesen 1999; Inouye 2020). The abundance of muscid flies declined by 80 percent over eighteen years at the Zackenburg field study sites in Greenland, an arctic site where muscid flies are dominant pollinators and one of the most diverse families of insects (Loboda et al. 2017). Such dramatic declines suggest that incorporating more fly taxa into long-term monitoring programs in northern regions should be prioritized, even though identification for many dipteran families is challenging. Long-term responses can vary among families, and order-level abundance trends may not be representative of individual taxa (Høye et al. 2021).

Among epigeic arthropods, ground beetles (Carabidae) are the most commonly used bioindicators for long-term monitoring programs like the National Ecological Observation Network (Hoekman et al. 2017), but rove beetles (Staphylinidae) represented the overwhelming majority of beetle abundance and diversity in our study area. Unfortunately, rove beetles have been eschewed by many studies despite their sensitivity to vegetation gradients and disturbances at multiple spatial scales (Klimaszewski et al. 2021). Other overlooked but potentially informative taxa include mites (Acari) and springtails (Collembola). New resources for identification, including improved taxonomic keys and DNA barcode libraries, may make these groups more accessible for future long-term monitoring programs. Sampling diverse taxa can lead to more robust conclusions than sampling any one taxon alone (Bergeron et al. 2013; Ernst, Loboda, and Buddle 2016; Høye et al. 2021).

Conclusions

We characterized arthropod communities along a montane elevation gradient to assess how these communities may respond to climate-driven shifts in habitat type boundaries. We showed that community composition clearly differs between habitat types, and these changes were strongly correlated with temperature and vegetation cover. If current warming trends continue, we predict that many forest- and shrub-associated arthropods will be able to track vegetation boundaries upslope, although idiosyncratic responses for some species have the potential to change community structure and function. The rapidly shrinking alpine tundra contains a high proportion of specialist arthropod species that are unlikely to thrive in a warmer, shrubbier environment, and upslope expansion by these species out of tundra and into unstable, rock, scree, and talus, with relatively few plants, is also unlikely. Recent reports of global insect declines have underscored the need for long-term data sets, which are rare at high latitudes where the effects of climate change are most pronounced. We established critical baseline data, established a transect system for long-term monitoring, and identified indicator taxa to help meet this goal.

Acknowledgments

We thank interns, students, and volunteers who assisted with field and laboratory work: Juan Diego Aguilar, Felix Bruner, Kelvin Chen, Carolene Coon, Kathryn Daly, Laura Larson, Renee Nowicki, Ashley Mocorro-Powell, Evelin Preciado, Amber Roberts, Alan Roos, Tazheem Rubio, Zachary Snelson, Rob Stewart-Rogers, Jayce Williamson, Mark Winterstein, Mary Wyatt, and Duncan Ziegler. We benefited from the expertise of many taxonomists for species identifications and are grateful to Adam Brunke (Coleoptera: Staphylinidae: Quedius), Shawn M. Clark (Coleoptera: Chrysomelidae), Joyce Cook (Coleoptera: Leiodidae), Hume Douglas (Coleoptera: Elateridae), Michael Geiser (Coleoptera: Cantharidae), Jason Gibbs (Hymenoptera: Lasioglossum), Terry Griswold (Hymenoptera: Andrena, Osmia), Paul Johnson (Coleoptera: Byrrhidae, Elateridae), Dave Kavanaugh (Coleoptera: Carabidae: Nebria), Jan Klimaszewski (Coleoptera: Staphylinidae: Aleocharinae), David Langor (Coleoptera: Carabidae), Ted MacRae (Coleoptera: Buprestidae), Kevin Moran (Diptera: Syrphidae), Alfred Newton (Coleoptera: Staphylinidae), Stewart Peck (Coleoptera: Leiodidae), Georges Pelletier (Coleoptera: Cryptophagidae), Darren Pollock (Coleoptera: Pythidae), Alexander Ryvkin (Coleoptera: Staphylinidae: Steninae), Alexey Shavrin (Coleoptera: Staphylinidae: Omaliinae), Joey Slowik (Araneae), Andrew Smith (Coleoptera: Scarabaeidae), and Margaret Thayer (Coleoptera: Staphylinidae: Omaliinae). We thank Scott Armbruster, Greg Breed, and Teresa Hollingsworth for their research guidance.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This work was supported by the U.S. National Park Service (Study #DENA-00925), The Alaska Entomological Society, and the Entomological Society of America SysEB Section.