Empirical landscape genetic comparison of single nucleotide polymorphisms and microsatellites in three arid‐zone mammals with high dispersal capacity

Abstract Landscape genetics is increasingly transitioning away from microsatellites, with single nucleotide polymorphisms (SNPs) providing increased resolution for detecting patterns of spatial‐genetic structure. This is particularly pertinent for research in arid‐zone mammals due to challenges associated with unique life history traits, such as boom‐bust population dynamics and long‐distance dispersal capacities. Here, we provide a case study comparing SNPs versus microsatellites for testing three explicit landscape genetic hypotheses (isolation‐by‐distance, isolation‐by‐barrier, and isolation‐by‐resistance) in a suite of small, arid‐zone mammals in the Pilbara region of Western Australia. Using clustering algorithms, Mantel tests, and linear mixed effects models, we compare functional connectivity between genetic marker types and across species, including one marsupial, Ningaui timealeyi, and two native rodents, Pseudomys chapmani and P. hermannsburgensis. SNPs resolved subtle genetic structuring not detected by microsatellites, particularly for N. timealeyi where two genetic clusters were identified. Furthermore, stronger signatures of isolation‐by‐distance and isolation‐by‐resistance were detected when using SNPs, and model selection based on SNPs tended to identify more complex resistance surfaces (i.e., composite surfaces of multiple environmental layers) in the best‐performing models. While we found limited evidence for physical barriers to dispersal across the Pilbara for all species, we found that topography, substrate, and soil moisture were the main environmental drivers shaping functional connectivity. Our study demonstrates that new analytical and genetic tools can provide novel ecological insights into arid landscapes, with potential application to conservation management through identifying dispersal corridors to mediate the impacts of ongoing habitat fragmentation in the region.


| INTRODUC TI ON
To conserve biodiversity, it is essential to preserve the evolutionary processes that support it, such as dispersal, mating, gene flow and selection (Latta, 2008). Incorporating dispersal knowledge into conservation planning is fundamental as this represents where and when species move in the landscape (Driscoll et al., 2014). However, deriving dispersal estimates through methods such as capturerecapture or telemetry is expensive and difficult in non-abundant species (Waits et al., 2015). Furthermore, the distinction between long-distance dispersal events that lead to gene exchange (realized dispersal), crucial for maintaining functional connectivity between populations, versus home-range movements is challenging to ascertain (Jordano, 2017). A valuable and important proxy for measuring functional connectivity is the empirical estimation of gene flow using genetic data since this is equivalent to measuring realized dispersal (Manel et al., 2003;Whitlock & McCauley, 1999).
Landscape genetics aims to determine how realized dispersal is influenced by the surrounding environment. It operates by examining genetic variation within heterogenous landscapes to explicitly quantify the effect of landscape composition and/or matrix quality (the space separating habitat patches) on an organism's dispersal, gene flow and/or spatial-genetic structure (Manel et al., 2003).
Microsatellite markers have previously been the major tool for such research (Storfer et al., 2010), capitalizing on short nucleotide motifs that are repeated in tandem a variable number of times. Due to high mutation rates, they possess a high information content per locus (Morin et al., 2004), although panels are often constrained to just tens of markers and represent a small fraction of the genome (O'Leary et al., 2018). Alternatively, single nucleotide polymorphisms (SNPs) are abundant and widespread throughout the genome (Morin et al., 2004). Though SNP markers are bi-allelic and give lower information content per locus, this can be offset by the large number generated (often thousands to tens of thousands) (Andrews et al., 2016).
With appropriate sampling, microsatellites have revealed patterns of functional connectivity and landscape barriers to dispersal across a variety of species and ecosystems (Emaresi et al., 2011;Munshi-South, 2012;Trénel et al., 2008). However, evidence suggests that SNP markers have higher accuracy and power to detect individual, population and species level patterns of genetic structure (Kim & Roe, 2021;Sunde et al., 2020). SNPs consistently outperform microsatellite markers in comparative studies analyzing population structure and assignment methods, specifically for finer-scale population genetic structure or species with high levels of gene flow (Jeffries et al., 2016;Puckett & Eggert, 2016;Viengkone et al., 2016).
However, to our knowledge, no studies have yet compared findings between marker types in relation to landscape genetic isolation-byresistance hypotheses (IBR; where dispersal is influenced by the degree of landscape resistance) (McRae, 2006).

Advances in genetic markers, technology and landscape genetic
methods provide opportunities to resolve patterns in cryptic landscapes, such as the topographically complex Pilbara region, situated in the Australian arid biome. The Pilbara is a biodiversity hotspot that supports rich faunal diversity including both endemic and widespread mammals ), yet the functional connectivity in Pilbara mammals is poorly resolved. This is due to several factors such as the boom-bust population dynamics that many aridzone mammals possess, as well as long-distance dispersal capacities to overcome sharp ecological gradients (Dickman et al., 1995;Kelly et al., 2013). The few genetic studies in the region are based on microsatellite and mitochondrial markers and reveal low genetic structure (Hohnen et al., 2016;Levy et al., 2019;Umbrello et al., 2020).
Threats including resource extraction, grazing pressure, and inappropriate fire regimes all impact habitat connectivity in the Pilbara (Cramer et al., 2016), highlighting the need to understand functional connectivity in the region.
Here, we assess spatial-genetic structure for three small grounddwelling mammals (body weight <15 g) adapted to arid environments, including a carnivorous dasyurid marsupial: Ningaui timealeyi, and two native rodents: the western pebble-mound mouse, Pseudomys chapmani and the sandy inland mouse, Pseudomys hermannsburgensis. While P. hermannsburgensis is widespread across most of arid Australia, both P. chapmani and N. timealeyi are Pilbara endemics. Both P. hermannsburgensis and N. timealeyi are habitat generalists, although P. hermannsburgensis shows a slight preference for sandy substrates (Gibson & McKenzie, 2009). Conversely, P. chapmani is a habitat specialist associated with rocky substrates, requiring small, uniform pebbles to construct mounds (Start et al., 2000).
A previous study found limited evidence that the Pilbara landscape influenced patterns of genetic structure in these species (Levy et al., 2019). With new, high-resolution genetic and spatial data, we investigate patterns of functional connectivity in the Pilbara using microsatellites and SNPs and qualitatively compare the ability of each data set to resolve (1) population genetic structure and potential physical barriers to dispersal (isolation-by-barrier; IBB); and (2) the role of dispersal capacity (isolation-by-distance; IBD) and specific landscape attributes (isolation-by-resistance; IBR: aridity, soil moisture, substrate, topography, distance to water, vegetation, and/ or fire) in facilitating or restricting realized dispersal. We explore how vast and dynamic arid landscapes shape the spatial-genetic structure of arid-zone species with high capacity for dispersal, and how new analytical and genetic tools can provide novel ecological insights for conservation.

| Study area
The Pilbara bioregion covers an extensive 179,000 km 2 and is di- consisting of extensive marsh and flood-out zones, and Roebourne is comprised predominantly of sandy coastal plains . The two main bioclimatic zones overlapping the Pilbara include semi-tropical and arid climates (Sudmeyer, 2016).

| Genetic data sets
Microsatellite data for N. timealeyi (Nt), P. chapmani (Pc) and P. hermannsburgensis (Ph) were obtained from Levy et al. (2019) (Ph). A total of 100-500 ng of genomic DNA was sent to DArT for library preparation and sequencing.
Library preparation by DArTseq™ follows a reduced representation method with enzyme digestion (here PstI and SphI) followed by sequencing on an Illumina Hiseq 2500 (Nt: medium density sequencing at 1.2 million reads; Pc and Ph: high density sequencing at 2.5 million reads) (Cruz et al., 2013;Kilian et al., 2012). Read assembly, quality control and SNP calling was carried out through DArT's proprietary software (Melville et al., 2017).
The raw SNP data sets (Nt = 36,899,Pc = 45,733 and Ph = 139,916 SNP loci) were filtered in R version 4.1.2 (R Core Team, 2021) using a custom R script (Shaw, 2022;Shaw et al., 2022) with functions from dartr (Gruber et al., 2018) and SNPRelate (Zheng et al., 2012) packages. Filtering thresholds were determined by visualizing the raw data (see Appendix S1) and filtering was completed sequentially in order of appearance in text. First, individual call rate filtering resulted in the removal of one N. timealeyi individual (call rate <0.5), and two P. hermannsburgensis individuals (call rate <0.7). Locus quality filters were then applied to each data set based on thresholds for missing data (Nt & Ph: 5%, Pc: 4%), average total read depth (Nt: between 20 and 100, Pc: between 15 and 25, Ph: between 20 and 40), repeatability average (95%), and minimum minor allele frequency (using thresholds that removed singletons: Nt & Ph: 0.025, Pc = 0.05). Multiple SNPs per sequence were removed followed by linkage disequilibrium pruning. We calculated pairwise relatedness (Wang, 2002) in the R package related (Pew et al., 2015), using a F I G U R E 1 Map of the focal study site located in the Pilbara region of north-west Western Australia. The greyscale Digital Elevation Model overlaid outlines regional topography. Nature reserves are indicated in purple; where NP = National Park and CP = Conservation Park, and the major Interim Biographic Regionalisation of Australia (IBRA) subregions of the Pilbara are highlighted in orange.
relatedness threshold of 0.24 to remove additional samples (Nt = 2, Pc = 14 and Ph = 4) to avoid biasing the genetic analyses with closely related individuals (i.e., half-siblings and above) (Wang, 2018). We included an additional filter to remove loci not in Hardy-Weinberg Equilibrium (HWE) for IBD and population genetics summary statistics using dartr to remove SNPs that significantly deviated from HWE assumptions, with a Bonferroni correction. HWE filtering was carried out within genetic clusters (i.e., populations) identified through IBB analysis, described below.

| Isolation-by-barrier: genetic clustering
The presence of barriers to dispersal in the three focal species was investigated using multiple analyses. Firstly, we ran a Principal Coordinate Analysis (PCoA) (Legendre & Legendre, 2012) in dartr to identify natural genetic clusters in the data (based on Euclidean genetic distances). As a preliminary step, we visualized temporal patterns to ensure genetic structure was not related to sample collection year (Appendix S4). Next, we visualized spatial patterns in the data, and these results helped guide the maximum value for K (the number of ancestral populations) when using the TESS3 algorithm in the R package tess3r (Caye et al., 2016). As opposed to Bayesian clustering programs like STRUCTURE that utilize Markov chain Monte Carlo (MCMC) methods (Pritchard et al., 2000), tess3r estimates individual ancestry coefficients based on sparse non-negative matrix factorization algorithms (sNMF) taking geographic information into account (Caye et al., 2016). This algorithm produces similar results to Bayesian clustering methods while being substantially faster (Frichot et al., 2014). Unlike STRUCTURE and related models, this approach does not rely on assumptions such as linkage equilibrium and HWE in ancestral populations (Caye et al., 2016). Given our large geographic and temporal spread of samples, we deemed this the most appropriate model for our study system.
We tested K values one through seven, with 50 repetitions for each value and the maximum number of optimisation iterations set to 200. We used the default settings for the remaining parameters and masked 10% of the data to use for the cross-validation. The best performing run with the lowest root mean squared error (RMSE) was presented, with the best value for K decided based on the presence of a plateau or change in slope in the cross-entropy criterion.

| Isolation-by-distance and summary statistics
We performed tests for IBD and calculated population genetic summary statistics for each genetic cluster identified in the IBB analysis for both microsatellite and SNP data sets. Individuals were assigned to a genetic cluster if the corresponding admixture coefficient proportion was ≥0.7, while individuals were excluded from IBD and summary statistics if they fell below this threshold. We chose this somewhat arbitrary threshold to enable the representation of discrete populations for comparison and these individuals were not excluded from any other analyses. We also removed three P. hermannsburgensis individuals located on islands from IBD analysis and summary statistics, given the low likelihood that mainland and island individuals are randomly mating. IBD was investigated with Mantel tests in GenAlEx 6.5 (Peakall & Smouse, 2006Smouse et al., 1986) using 999 permutations to determine the significance of the Mantel correlation coefficient (Rxy). Summary statistics were then calculated in GenAlEx 6.5 within each genetic cluster, including the number of alleles (N a ), observed heterozygosity (H o ), expected heterozygosity (H e ) and the fixation index (F). When more than one genetic cluster was detected, we also calculated F IS and F ST .

| Isolation-by-resistance
We collated and derived high resolution spatial layers to test IBR landscape genetic hypotheses (Table 1; Appendix S2). We sought to represent aspects of aridity (aridity indices-ADI, ADX, and ADM), landscape productivity (relative soil moisture indices-SOMO29 to 30), substrate (clay, sand, silt, and coarse fragments-CF), topographic features (weathering intensity index-WII, vector ruggedness measure-VRM, digital elevation model-DEM), watercourses (Euclidean distance to water-WAT), vegetation (spinifex density index-SPIN, persistent forest cover-FOR), and fire (fire frequency-FF). The justification for selecting these variables along with definitions, data sources, and relevant citations can be found in Table 1. All layers were aggregated to a 5 km 2 resolution due to the long-range dispersal capacities of the focal species. For example, both Pseudomys species have displayed long-range movements of several kilometers, with P. chapmani individuals found up to 2 km from neighboring mounds and P. hermannsburgensis recorded moving up to 14 km in a 2-week period (Dickman et al., 1995;Start et al., 2000). Evidence from dasyurid species suggests that longrange dispersal of several hundreds of meters to several kilometers is also likely for N. timealeyi (Dickman et al., 1995).
The parameterisation of resistance surfaces within landscape genetic analyses has traditionally relied on subjective "expert opinion" which can sometimes lead to inaccuracy (Liu et al., 2018).
Furthermore, researchers generally assume a linear relationship between continuous variables and genetic distance despite this often not being the case . For these reasons, we used a genetic algorithm to parameterise resistance surfaces (i.e., relationship between pairwise genetic and effective distances) and maximum value (i.e., cost ratio) through optimizing for the best transformation with no a priori assumptions (Peterman et al., 2014), and fit this relationship using linear mixed effects models (described below). This was implemented in the R package ResistanceGA (Peterman, 2018).
Given these analyses are sensitive to contemporary patterns of gene flow, we removed four samples (Nt = 1, Ph = 3) with missing date information from this analysis. We also removed individuals located on islands (Ph = 3) from the analysis, as our aim was to explore how terrestrial landscape attributes influence gene flow (rather than determining whether the ocean is/is not a barrier to dispersal). Shirk TA B L E 1 Spatial layers used to test isolation-by-resistance, with proposed mechanisms and hypotheses detailed in the justification column (see Appendix S2 for layer visualization and further information).

Justification for hypotheses
Included in final variable set
Areas of higher aridity inhibit primary vegetation growth and invertebrate communities (Walsberg, 2000). We predict this will limit gene flow through reduced protection from predators and limiting food resources while dispersing.
Areas of higher soil moisture stimulate primary vegetation productivity and resource availability (Berndtsson et al., 1996). We predict this will facilitate gene flow by providing protection from predators and food resources while dispersing.
Clay is an important predictor of Nt occurrence; Ph prefers increasingly sandy substrates; Pc prefers rocky substrates (Ford & Johnson, 2007;Gibson & McKenzie, 2009). We predict that species' preferred substrate will be easier to move through, thus facilitating gene flow.

MSAT:
Nt WII, VRM, DEM WII: Weathering Intensity Index -describes regolith properties (low = unweathered outcrops, high = areas with clays and sands); VRM: Vector Ruggedness Measure -measure for terrain complexity (independent of slope); DEM: digital elevation model. Nt associated with rugged topography; Ph commonly occupies gentle topography; Pc displays preference for eroded, hilly areas of unweathered bedrock (Ford & Johnson, 2007;Gibson, 2011). We predict that species' preferred topography will be easier to move through, thus facilitating gene flow.
Other mammals in the Pilbara have shown positive association with permanent water, likely due to higher quality habitat (Moore et al., 2019). We predict proximity to water will facilitate gene flow by providing food resources and protection from predators. Higher vegetation density (e.g., Triodia) provides protection from predators for small mammals (Moseby et al., 2016).
We predict this will facilitate gene flow by providing increased protection from predators while dispersing.  (Goslee & Urban, 2007). Effective distances were calculated within the ResistanceGA package based on random-walk commute times, equivalent to CIRCUITSCAPE "resistance distance" (Klein & Randić, 1993;McRae et al., 2008).
The creation of optimized composite resistance surfaces through ResistanceGA can be applied to both categorical or continuous rasters and can be performed either independently or simultaneously across all raster layers. This is achieved by fitting linear mixed effects models with a maximum likelihood population effects parameterisation (MLPE) (Clarke et al., 2002) to the pairwise genetic data (response) and effective distances (predictor). During optimisation, models are compared based on an objective function (we used the default option, log-likelihood) across different transformations and parameters over "generations" until there is no improvement, thus indicating the best optimized surface.
Before calculating the multi-surface optimisations, we reduced collinearity in the raster data by removing correlated variables (Spearman's |r s | > .7; Appendix S3), by running single surface optimisations and selecting the top ranked surface in correlated sets (according to Akaike's Information Criterion corrected for small sample size; AICc) (Akaike, 1974). We included IBD (where each pixel in the resistance surface is given a value of 1) and a null model (intercept only) in the model selection, which is built into the ResistanceGA package. Surfaces that ranked lower than, or within 2 AICc of the IBD or null models in the single surface model selection described above were excluded from multi-surface optimisation, as they performed no better at modeling functional connectivity than the alternate or null hypotheses (Burnham & Anderson, 2002). Next, we performed multi-surface optimisation using the "all_comb" function in ResistanceGA on a maximum of four combined surfaces. We conducted 1000 bootstrap iterations across random subsets of 75% of the total data to calculate the percentage of iterations where surfaces were ranked as the top model (similar to model weight) (Burnham & Anderson, 2002). This provides an indication of the level of support for each surface and whether particular sample subsets are disproportionately influencing model results (e.g., samples collected at different time points). Model performance was assessed by visualizing residuals based on a simulation approach using the DHARMa R package (Hartig, 2022). Burnt areas provide less protection from predators (Moseby et al., 2016). We predict areas that are frequently burnt will limit gene flow through increased predation while dispersing.

| Isolation-by-barrier: genetic clustering
We found no indication of temporal genetic structuring in our data sets (Appendix S4). There was no evidence for IBB within the rodent species using both marker types, while for N. timealeyi, SNPs resolved two genetic clusters that were not detected when using microsatellite markers (Figures 2 and 3). Although this sample was also collected at the earliest timepoint (1988), we took this to represent genetic differentiation between island and mainland individuals given the lack of temporal clustering in the rest of the data set (Appendix S4). As with N. timealeyi, these same patterns were not resolved when using microsatellites ( Figure 2). Furthermore, tess3r analysis did not detect population genetic structure at either SNPs or microsatellite markers, suggesting the Pilbara represents one genetic cluster (K = 1) for both Pseudomys species (Figure 3a).

| Isolation-by-distance and summary statistics
Mantel tests and summary statistics were calculated across the total data set, as well as within the two genetic clusters (or "populations"; K = 2) for the N. timealeyi SNP data set identified in IBB analyses, excluding the 23 admixed individuals that did not assign to either population. For the Pseudomys SNP data sets, and the microsatellite data sets for all three species, analyses included all individuals as a single genetic population (K = 1).
We detected significant IBD (p < .05) within almost all populations identified in IBB analysis for each marker type for all three species (except for N. timealeyi cluster 1; Figure 4). However, while the weak positive relationship between genetic and geographic distance was consistent between marker types for P. hermannsburgensis

| Isolation-by-resistance
The number of variables retained following single-surface optimisation varied from one to eight (Table 1), with microsatellite-based model selection consistently resulting in fewer variables being retained than for SNPs (i.e., variables performed worse than or equivalent to IBD or the null model). For the remaining variables, multi-surface optimisation generated single and composite resistance surfaces that performed better in the model selection than the null model and IBD, for all three species across both marker types.
However, SNP data sets better differentiated between the best models (models ≤2 ΔAICc from the top-ranked model) and both IBD and null models, with ΔAICc for SNP IBD and null models 8 to 58 times greater (i.e., further from the top model) than for those using microsatellites (   ΔAICc for the second-best model = 6.611; Table 3; Appendices S5-S7), or clay (weight = 0.798, AICc = 4698.857, marginal R 2 = .025;

F I G U R E 4
Mantel plots for each genetic cluster (or "population") across marker types and species, where Rxy is the Mantel correlation coefficient and p indicates the significance using 999 permutations. Note that although only one genetic cluster was detected for N. timealeyi using the microsatellite data set, here we present the results for two clusters to match those identified using SNPs (for comparison), as well as for the total data set (i.e., where admixed individuals have not been removed). ΔAICc for the second-best model = 4.582; Table 3; Appendices S5-S7), respectively. The best supported SNP model ranked first in 92.6% of bootstrap iterations, compared to 79.3% for the microsatellite model (Table 3; Appendix S8). In both cases, landscape resistance increased as either ruggedness or clay content increased (Appendix S9). Although the SNP and microsatellite genetic TA B L E 3 MLPE model selection across all species and marker types, for models performing within 2 ΔAICc of the top-ranked model, as well as the isolation-by-distance (IBD) and null models for comparison (full model summaries and diagnostic plots can be found in Appendices S5-S7).

F I G U R E 5
Forest plot displaying MLPE model estimates, 95% confidence intervals and significance (asterisks) for the three different species. Models were fit with optimised single/composite resistance surfaces as predictor variables and Euclidean genetic distance as the response variable and built using microsatellite (orange circles) versus SNP datasets (blue triangles).

Microsatellites
SNPs response data selected different surfaces, both found that nonsandy substrates (rocky/rugged terrain and clay) increased landscape resistance, with the magnitude of this effect approximately 5.5 times greater when using the SNP dataset ( Figure 5).

| DISCUSS ION
Increasingly the field of conservation genetics is transitioning to the use of genomic data to understand patterns of connectivity in wild populations. While SNP markers appear to provide increased resolution for detecting population genetic structure, so far there has been a lack of studies to investigate the performance of SNPs relative to microsatellites for empirical landscape genetic analyses. Here we provide a case study using these two marker types to evaluate three explicit landscape genetic hypotheses (IBD, IBB and IBR) in a suite of small arid-zone mammals possessing high dispersal capacities. In general, SNP markers provided additional resolution in detecting subtle genetic structuring in IBB analyses, particularly for the dasyurid, and resulted in stronger patterns of IBD and IBR in both rodent and dasyurid species. However, the large dispersal capacity and the dynamic nature of arid landscapes present specific challenges for IBR analyses, which we discuss below.

| Isolation-by-barrier: SNPs versus microsatellites
Here, using SNP data, we identified genetic structure in the dasyurid species, N. timealeyi, that was not detected with microsatellite data using the same individuals (in this study and in Levy et al., 2019).
Our results add to a growing body of work suggesting SNPs provide higher resolution for population genetic analyses compared to microsatellites, particularly in species with weak population structure (Jeffries et al., 2016;Puckett & Eggert, 2016;Viengkone et al., 2016). Similar to our results, Camacho-Sanchez et al. (2020) compared both marker types in two amphibians, concluding that SNP data sets with large numbers of loci are more reliable at identifying population genetic structure at large spatial scales (~500,000 square kilometers). Furthermore, a synthesis of studies using both marker types showed that SNPs were either equivalent or outperformed microsatellites at detecting population genetic structure, suggesting this pattern is broadly representative, rather than study or context specific (Sunde et al., 2020).
Single nucleotide polymorphism markers indicated two genetically distinct groups in N. timealeyi, although given the large degree of admixture, these clusters may be better described as represent- Alternatively, given that SNPs have a slower mutation rate than microsatellites (Morin et al., 2004), this pattern may also reflect the ancestral genetic signature present prior to aridification in the midlate Pleistocene (approximately 15-25 kya). This period saw rivers transition from perennial to ephemeral flows and the expansion of drought resistant flora (Byrne et al., 2008). Thus, a third scenario of reconnection of previously separated refugial populations may also explain the substantial admixture and low differentiation between these genetic clusters. In fact, Umbrello et al. (2020) found evidence of population expansion in six small dasyurids across the Pilbara since the mid-late Pleistocene and the Last Glacial Maximum (LGM) and proposed that this followed the increased availability of arid habitat. Refugial separation prior to population expansion after the LGM has also been detected in sea spiders (Soler-Membrives et al., 2017), mussels (Cunha et al., 2011) and ants (Xun et al., 2016), with weak differentiation reflecting the loss of refugial genetic structure over time due to high dispersal capacities.
In contrast, we were not able to detect evidence for population genetic structuring within the Pilbara landscape for the two native Pseudomys species across both SNP and microsatellite data sets, suggesting a lack of landscape barriers to dispersal. However, SNPs were still able to resolve some subtle patterns not detected with microsatellites (e.g., PCoA groupings of the P. hermannsburgensis island individual and the subtle western and north-south groupings for P. chapmani). The weak clusters detected in the SNP P. chapmani data set may reflect the accumulation of positive spatial-genetic structure driven by the sociality of the species (i.e., family groups within pebble mounds) Ford & Johnson, 2007).
However, this structure was too weak to be detected with our sampling strategy (spatially dispersed individuals) and the clustering analysis. Several other genetic studies in rodents also find low population structure even in the presence of major landscape barriers or considerable landscape heterogeneity (Gauffre et al., 2008;Vega et al., 2007). This is likely because irruptive boom-bust population dynamics obscure any signals of population structure.

| Isolation-by-distance and isolation-byresistance: SNPs versus microsatellites
Few studies have evaluated the ability of microsatellites versus SNPs to detect IBD and we are not aware of any that have compared results between marker types for identifying IBR. In a comparative study using RADseq SNPs and microsatellites, Jeffries et al. (2016) identified a stronger signature of IBD from the SNP data than from microsatellite data sets, suggesting this may be due to the mutational processes of the markers. Similarly, we detected significant IBD for all species and marker types, and this signature was stronger when using SNPs in some cases. For example, while both marker types showed weak IBD for P. hermannsburgensis, the magnitude of the correlation increased when using SNPs for N. timealeyi and P. chapmani. Perhaps the increased power provided by more loci, coupled with the slower mutation rate of SNPs was able to resolve this subtle pattern, suggesting that these species have more constrained dispersal capacities than P. hermannsburgensis. For IBR, we found that SNPs tended to resolve more complex resistance surfaces (i.e., composite surfaces of multiple environmental layers) than microsatellites, potentially reflecting the increased power of large SNP panels to detect subtle and complex patterns of functional connectivity. SNP models also revealed a stronger effect of landscape resistance on genetic distance and tended to better differentiate between the top models and the alternate IBD hypothesis, adding to the body of evidence arguing that SNPs provide better resolution for questions that require individual-level genetic information, such as relatedness, individual identification and fine-scale genetic structure (Sunde et al., 2020).
High resolution spatial data in combination with sophisticated landscape resistance modeling revealed additional detail on the landscape characteristics influencing functional connectivity in our target species than detected in Levy et al. (2019), and the identified drivers of connectivity were largely consistent between marker types. When directly compared via model selection, we found that IBR hypotheses outcompeted the alternate hypothesis of IBD in all cases. However, these results should be interpreted carefully since the high power inherent to large pairwise data sets, combined with correlations between competing IBD and IBR models, can result in low model selection accuracy (Shirk et al., 2018). The large dispersal capacity and boom-bust dynamics of arid-zone mammals makes this issue particularly pertinent for our study, given that functional connectivity would likely be approaching IBD for many species.
Given these issues, we attempted to reduce model selection error by following best practice recommendations, including using linear mixed-effects models fit with MLPE and by transforming resistance surfaces to satisfy assumptions of linearity, as this approach has been shown to outperform other regression methods (Shirk et al., 2018). We also used individual-genetic distance, which is sensitive to contemporary genetic structure (Shirk et al., 2017), since population-level analyses are less representative of species that are continuously distributed. Finally, we used a data-driven approach to parameterising resistance surfaces based on spatial layers that are biologically plausible, to tease apart competing hypotheses and determine the most likely characteristics contributing to landscape resistance.
Simulations can be used to undertake power analyses and evaluate findings (i.e., to determine whether marker panels have the power to detect genetic patterns in specific systems or scenarios).
Simulation tools such as CDPOP (Landguth & Cushman, 2010) and HexSim (Schumaker & Brookes, 2018) have contributed greatly to this goal, however, it can be difficult to parameterise simulations, particularly in relatively understudied systems such as the arid landscape presented here, due to complex and unknown species' demography.
In particular, it is not yet feasible to simulate landscapes where genetic structure plays out over such a vast scale and there is still much work to be done to develop landscape genetic tools to help us understand the interplay between boom-bust dynamics, temporally and spatially dynamic refuges (common features of arid landscapes), and dispersal.
Further research on arid systems can provide greater mechanistic understanding of these patterns and processes. Although we cannot make a quantitative and comparative power analysis without further simulation testing in our study, we can interpret our results alongside the literature to determine biological plausibility and provide useful ecological insights for conservation management.
Using this approach across both marker types, we found that increasing ruggedness or rockiness facilitated landscape connectivity for N. timealeyi. This species is a habitat generalist, weakly associated with clay substrates and also found on rugged terrain (Gibson & McKenzie, 2009) gensis both showed top models that were orders of magnitude higher than the IBD (and null) models, this difference was slightly less pronounced for P. chapmani (particularly using the microsatellite data set), suggesting that functional connectivity is approaching IBD for this species. Still, our results indicated that increased rocky outcrops (lower weathering-SNPs) or increased rocky substrate (coarse fragments-microsatellites), and increased soil moisture (SNPs) facilitated landscape connectivity. This is biologically plausible since this species is a rocky habitat specialist, although our results suggest that due to the fragmented nature of rocky habitat in the Pilbara (Ford & Johnson, 2007), dispersal must often occur between these patches, a pattern also seen in other species found in rocky habitat in the Pilbara .

| CON CLUS IONS
There has been a rapid shift from microsatellite markers to SNPs in the fields of conservation and population genetics, although studies suggest that the major benefit of SNPs is not inherently about the marker type, but the number used (Sunde et al., 2020). We found that SNPs resolved subtle genetic structuring not detected by microsatellites, stronger signatures of isolation-by-distance and isolation-byresistance, and identified more complex resistance surfaces. While patterns of genetic structure were subtle, our study demonstrates that the use of SNPs, coupled with novel landscape genetics analyses, can provide new ecological insights into arid landscapes, although microsatellites were still able to identify similar (albeit more simplified) results.
Understanding subtle resistance patterns in highly permeable landscapes is not of obvious conservation concern (Shirk et al., 2018).
However, the Pilbara is substantially impacted by competing land uses, including mining and pastoralism (and the cumulative impacts of habitat clearance and fragmentation), as well as climatic cycles that drive dynamic drought and fire regimes (Cramer et al., 2016(Cramer et al., , 2022McKenzie et al., 2009). Even in large, panmictic populations, maintaining functional connectivity is crucial, as ongoing fragmentation can erode meta-population health (Umbrello et al., 2022). Thus, knowledge of factors driving connectivity is crucial for supporting resilient populations of both threatened and non-threatened species in Australia, and globally. data curation (lead); formal analysis (lead); investigation (supporting); methodology (lead); supervision (lead); visualization (equal); writing -original draft (supporting); writing -review and editing (lead).

ACK N OWLED G M ENTS
This study was part of an Australian Research Council Linkage Project Government and the Government of Western Australia. We gratefully acknowledge the many Traditional Owners of the Pilbara region; on whose country these data were collected. We are also grateful for data and support provided by Murdoch University and DBCA. In particular, we acknowledge the personnel involved with the Pilbara biological survey and Janine Kinloch for providing spatial data. Open access publishing facilitated by Australian National University, as part of the Wiley -Australian National University agreement via the Council of Australian University Librarians.

B EN EFIT-S H A R I N G S TATEM ENT
All data and R code are available from Github (https://github.com/ Robyn Sh/LandG en_AridM ammals) and have been archived in a Zenodo Digital Repository (Shaw, 2023).
This research is part of a collaborative project across academic, government and industry partners and addresses a priority concern by testing the methods being applied to the conservation of small mammals in the broader project. Benefits from this research accrue from the sharing of our data and results on public databases as described above.