Introduction
Mountainous tropical regions are hyperdiverse (Fjeldså et al., 2012; Rahbek, et al., 2019a; Rahbek, et al., 2019b). Besides climate heterogeneity producing different habitats in a patchy distribution, it has been argued that two additional factors particularly increase diversity in tropical mountains. First, because tropical taxa have narrow physiological tolerances to temperature, tropical mountain passes are more effective barriers to dispersal than those in temperate regions, thus promoting speciation through physical disruption of gene flow (Janzen, 1967; Polato et al., 2018; Sheldon et al., 2018). Second, coupling altitudinal gradients with tropical latitudes allows populations within tropical mountains to persist relatively in situ despite global climate fluctuations (Mastretta-Yanes, et al., 2018; Rahbek, et al., 2019a; Rahbek, et al., 2019b). These processes of isolation and long-term persistence have been widely used to explain diversification among tropical mountain peaks across the world (e.g., Fjeldså et al., 2012; He et al., 2019; Knowles, 2001; Mastretta-Yanes et al., 2018; McCormack et al., 2009; Uscanga et al. in review). However, recent evidence suggests that neutral processes within a single mountain could also play an important role in generating endemism (Bray & Bocak, 2016), but diversification at the single mountain scale and short geographic distances has seldom been explored.
Studying evolutionary processes in mountain ranges across fine spatial scales is challenging due to their complexity with regard to topography and geological and climatic history. However, the insular nature of sky-islands reduces the complexity of mountainous regions, thus providing more simplified systems to analyze evolutionary processes. Within true oceanic islands, habitat discontinuity has been demonstrated to have an important role in driving geographical diversification (García‐Olivares et al., 2019; Goodman et al., 2012). Additionally, Salces‐Castellano et al. (2020) have demonstrated that, when dispersal ability and climate tolerance are restricted, strong geographic isolation over distances of only a few kilometers can be found for multiple co-occurring arthropod species of an oceanic island. These results suggest that besides allopatric diversification between different islands, island systems could also promote high levels of intra-island geographical diversification, acting as a local-scale diversity source. If geographical diversification were occurring within individual sky-islands in this way, following the expectations of the neutral theory of biodiversity, it would be expected to find similar spatial patterns of differentiation from the haplotype to the community levels (Baselga et al., 2013).
Arthropod communities are ideal systems to test evolutionary processes at fine spatial scales within sky-islands because they are locally abundant and diverse and relatively easy to sample massively. They also harbour a broad diversity of groups with different dispersal abilities, e.g. including winged and non-winged species and varying body-sizes. However, taxonomic identification to the species level of rich communities of arthropods with traditional methods is challenging (Favreau et al,. 2006, Yu et al., 2012). In this sense, new high throughput sequencing approaches applied to the study of arthropods are revolutionizing the understanding of complex arthropod communities (Andujar et al., 2015; Arribas et al., 2016; Ji et al., 2013; Yu et al., 2012). Specially promising is the use of whole community metabarcoding (cMBC) for the bulk sequencing of the mitochondrial COI gene of mixed communities (Andújar et al., 2018). Recent improvements for denoising metabarcoding datasets (e.g., Edgar, 2016; Callahan et al., 2016) and evaluating the prevalence of sequencing errors and co-amplified pseudogenes (Andújar et al., 2020) raised the prospect of read-based, haplotype-level analyses with mitochondrial COI cMBC data, which represents a step change for the study of diversity patterns through whole-community genetic analyses (Andújar et al., 2018; Arribas et al., 2020).
Haplotype data from hyperdiverse arthropod communities can be used directly for analyses of genetic diversity, or aggregated into species-level entities for analyses of species diversity, allowing for the joint analysis of turnover (beta diversity) at multiple hierarchical levels. Local assemblages may diverge simply due to the lack of population movement which, when assessed for entire communities, results in a largely regular decay of community similarity with spatial distance for the typically neutral haplotype variation of the mitochondrial COI gene (Baselga et al., 2013). Under a scenario where dispersal constraints are a dominant driver of spatial variation in community structure, assemblage turnover at the species level should mirror these haplotype patterns, albeit at a higher level of similarity (Baselga et al., 2013; Baselga, Gómez-Rodríguez, & Vogler, 2015). This analytical approach has been exploited to determine whether the composition across multiple beetle taxa assemblages is predominantly driven by dispersal (Baselga et al., 2013; Baselga, Gómez-Rodríguez, & Vogler, 2015) and has been proposed as a useful way to compare relative dispersal constraints between lineages from different taxonomic groups (Gómez-Rodríguez et al., 2019, Múrria et al., 2017). These, and other studies using a multi-hierarchical approach, have focused from regional to continental-scale distances. However recent work has also exploited this framework to analyze community assemblage at finest geographic scales (<15 km) using cMBC data, while also allowing to explore hyperdiverse communities, like soil mesofauna (Arribas et al., 2020). Therefore, applying the cMBC approach to the hyperdiverse arthropod faunas of tropical mountains offers much potential for understanding their community structure and the processes that have shaped it.
Here, we evaluate fine-spatial community patterns within the arthropod fauna of a tropical sky-island forest, to better understand the roles of dispersal limitation and landscape features as drivers of the diversity found in tropical mountains. To do this, we performed a systematic sampling consisting of 840 pitfall traps across Abies religiosaforests within the Nevado de Toluca, a sky-island from the Transmexican Volcanic Belt (TMVB). We generated haplotype-level metabarcoding data for 42 arthropod communities distributed in sampling blocks separated from 50 m to 19 km, and evaluate patterns of richness, turnover and distance decay in community similarity at multiple hierarchical levels (haplotype, putative species (OTU) and supra-specific levels). As we are interested in the effect of ecological and topographic features on dispersal limitation, in addition to testing for the effect of isolation-by-distance (IBD), we also performed an isolation-by-resistance (IBR) analysis. This allows cost dispersal given by landscape features to be incorporated, yielding biologically more informative distance decay relationships than Euclidean distance alone (McRae, 2006). Our analytical framework thus allows us to assess the role of dispersal constraint within a local spatial setting on tropical mountain diversity.