Our knowledge of community assembly dynamics under multiple stressors is limited because of opposite niche-based processes, i.e., limiting similarity and habitat filtering could simultaneously occur, masking the overall patterns. Alpine biomes provide the ideal systems to explore the influences of co-occurrence processes as these communities usually face multiple stresses such as resources limitation and habitat constraints. However, the assembly processes of mammals in alpine areas have hardly been exclusively studied. Here, we aimed to address how different processes structured small mammal communities at the tree line transition zone, which represents one of the most distinct vegetation boundaries separating alpine from montane habitats. We compiled a regional dataset including species list, phylogeny, and functional traits from field collections across 18 mountains of southwest China and complemented them with published data sources. The traits were decomposed into different niche components to determine the respective effects of specific stressors. Phylogenetic and functional diversity indices representing evolutionary history, trait space, and pairwise species distance were calculated and compared with null expectations. Linear mixed-effect models were constructed to assess the increasing or decreasing tendencies of diversity values against increasing elevation. The results showed that phylogenetic and functional richness were strongly correlated with species richness, unlike the distance-based indices which were uncorrelated with species richness. There was no evidence found to support non-random phylogenetic or overall trait patterns. However, the resource acquisition niche tended to be more overdispersed (positive slopes), while the habitat affinity niche tended to be more clustered (negative slopes) as habitats became less productive and less vegetated. We conclude that limiting similarity and habitat filtering simultaneously structure small mammal communities in alpine areas. Altogether, the present study provides vital insights into the complexity of co-occurring assembly processes by niche decomposition, and highlights the importance of considering different diversity dimensions when assessing community structure.

2.1 Study area and field sampling

The study region is located in the southwest mountains of China, a global biodiversity hotspot (Mittermeier et al. 2011), and covers the Three Parallel Rivers of Yunnan Protected Areas (Fig 1). Four north-south trending mountain ranges dominate the region's topography (from west to east, Gaoligong, Nushan, Yunling, and Shaluli). Eighteen prominent mountains (or sky-islands) were selected as study sites to adequately represent the study region's extent (Appendix 1). The tree line in the region represents one of the most distinct vegetation boundaries separating alpine from montane habitats (Körner 2013; Testolin et al. 2020). Across the study region, the alpine tree line is located between ca. 3 800 – 4200 meters above mean sea level (m a.s.l.), and generally rises from south to north but varies between sites due to several abiotic or anthropogenic influences (Wang et al. 2013). The alpine habitats in the region experience a moderate decline in annual precipitation from west to east (Sherman et al. 2008), with a 6.7 ℃ mean annual temperature at the tree line (Wang et al. 2013). The vegetation from below to above the tree line rapidly changes from dark coniferous forests to dwarf shrubs, meadows, and scree (Fig. 3; Sherman et al. 2008).

The small mammal sampling was conducted from 2013 to 2018 following a standardized field method (See Song et al. 2020 for details). Basically, three elevational transects (from 200 m below to 200 m above the tree line at 200 m intervals) were established on each mountain, representing the environmental transition from montane to alpine habitats. We chose this elevation range because it covered most small mammal species' strongest habitat transition zone before their upper distribution limits. The highest transect was missing in seven mountains because their summit merely transcended the tree line, and the remaining elevations beyond the tree line were not enough for another transect (＜ 200 m). For each captured small mammal, measurements of body weight (BW), head-body length (HB), tail length (TL), hindfoot length (HF), and ear length (EL) were recorded in the field. From each site (mountain), at least one stuffed voucher and skull specimen was prepared for each morphologically discriminable group. Liver and muscle tissues were collected from each specimen, dehydrated, and preserved in 2ml vials (81-8204, Biologix Inc., Shandong, China) with 99.7% ethanol for DNA extraction. All samples were identified by morphological comparison and DNA barcoding based on the mitochondrial Cytochrome b (CYTB) gene (Borisenko et al. 2008). All collections are deposited at Kunming Institute of Zoology, Chinese Academy of Sciences.

2.2 Reginal species pool

The regional species pool crucially influences conclusions on local communities’ assembly patterns and processes (Weiher and Keddy 1995). Species from the regional pool should have the potential to colonize, establish, or inhabit the focal communities depending on different ecological processes (Lessard et al. 2012) and be adequately detected in the field. The regional small mammal species pool was primarily derived from Wen et al. (2016). This data set contains county-level species occurrence records covering our study sites. To avoid eliminating influences of potential habitat filtering in alpine communities a priori (Lessard et al. 2012), we extracted species information from the counties where the 18 studied sites were located (Fig. 1), aiming to represent all small mammal species along the entire elevational gradient (from the lowest to highest elevations) within the region. Further, species that were unlikely to be captured by the standardized field sampling, such as marmot, flying squirrel, and hare, were removed. The final species pool included four Orders (Eulipotyphla, Lagomorpha, Rodentia, and Scandentia). Due to outdated taxonomical records for some species, we carefully checked the species list and replaced debatable names with the most recent taxonomic evaluations. New species records from the present study were noted (Song et al. 2021), and those that did not match any known species were assigned temporary taxonomic identifications. Species compositions in sites were Hellinger-transformed to reduce the effect of zero abundances (Legendre and Gallagher 2001).

2.3 Phylogenetic reconstruction

The CYTB gene has been found to be more accurate at reconstructing mammalian species' phylogenetic associations (Tobe et al. 2010) and is frequently used for small mammals’ phylogenetic studies (e.g. Bannikova et al. 2018; Koju et al. 2017; Liu et al. 2018). Total DNA was extracted from newly collected samples, and the CYTB gene was amplified and sequenced to rebuild the phylogeny. The DNA was extracted from 1-2 individuals of each morphologically distinct species from each site using the SQ Tissue DNA Kit (D6032, Omega Bio-tek Inc., Norcross, USA) following the manufacturer’s protocol. The CYTB gene was amplified and sequenced following Koju et al. (2017). The PCR primers for the genus Ochotona were L14724 and H15913 (Koju et al. 2017) and L14725_hsw1 (ATG ACA TGA AAA ATC ATC GTT GT) and H15915_hsw1 (TCY CCA TTT CTG GTT TAC AAG ACC) for other species. The sequenced CYTB products were assembled in Geneious R11.1 (https://www.geneious.com/) then aligned in MEGA X (Kumar et al. 2018) using MUSCLE (Edgar 2004). The resulting alignment was verified in BLAST (https://blast.ncbi.nlm.nih.gov/Blast.cgi) to handle sequencing errors and misidentifications. The alignment for further phylogenetic analysis comprised of one complete CYTB sequence of each newly captured species and the CYTB sequences of species not captured but downloaded from GenBank. The field numbers and GenBank accession numbers are provided in Appendix 2. The phylogenetic tree of the total species pool was reconstructed in the IQ-TREE web server (http://iqtree.cibiv.univie.ac.at/; Trifinopoulos et al. 2016) using Maximum Likelihood (ML) with 10 000 bootstrap replicates.

2.4 Traits used to quantify functional diversity

Traits were classified into morphology, resource acquisition, and habitat affinity combinations. Morphological characteristics, including body size and extremities (e.g., ears, limbs, and tail), are important to cold adaption (Alhajeri et al. 2020; Blackburn and Hawkins 2004), a prominent constraint in alpine areas. The morphological traits, represented by BW, HB, TL, HF, and EL, were measured in the field as numerical variables; BW and HB were represented by the mean values of recorded individuals, while TL, HF, and EL were transformed as proportions of HB (Du et al. 2017). Resource acquisition traits included activity cycle, trophic level, diet breadth, and diet composition extracted from the COMBINE dataset (Soria et al. 2021). Habitat affinity traits were represented by habitat breadth, habitat occupations, fossoriality, and foraging stratum obtained from Ding et al. (2022). Habitat breadth determines the number of available habitat types, and habitat occupation refers to the availability of specific habitat types for the species. Fossoriality and foraging stratum represented spatial separation in a given habitat, i.e., underground, above ground, or arboreal. Among all traits, diet composition and habitat occupation comprised 3 (vertebrate, invertebrate, and plants) and 6 (forest, shrublands, grassland, wetlands, rocky areas, and artificial) subdivisions, respectively. Diet composition included the percentage proportion of different food types (Soria et al. 2021), and habitat occupation included binary levels of availability (available =1 and unavailable = 0) of a particular habitat type for the species (Ding et al. 2022). Traits were classified into nominal, circular, or quantitative types (Magneville et al. 2021).

All statistical analyses and result illustrations were processed in the R 4.1.0 environment.

Initial data and R code were provided.

The uploaded files include:

 1. overall_trait.R: R code to calculate observed phylogenetic and overall trait functional diversity indices, null models, and plots.

2. morpho_trait.R: R code to calculate observed morphological trait functional diversity indices, null models, and plots.

3. die_trait.R: R code to calculate observed resource acquisition trait functional diversity indices, null models, and plots.

4. habit_trait.R: R code to calculate observed habitat affinity trait functional diversity indices, null models, and plots.

5. lme.R: R code to calculate correlation coefficients between diversity metrics and generating linear mixed-effects models.

6. pool.sp.csv: species by site matrix of the regional species pool

7. pooltree: netwick file of the phylogenetic tree used to calculate phylogenetic diversity and tests phylogenetic signals.

8. reduced_traits.csv: trait data of the regional species pool used to calculate functional diversity metrics.

9. tr_cat.csv: sheet used to classify the traits into nominal, circular, or quantitative types, which was required in the ‘funct.dist’ and ‘quality.fspaces’ in ‘mFD’ package. 

10. meta-redueced.csv: elevational gradients used to generate linear mixed-effects models.

11. README_file.txt: description for the variable in each data set.