Sex-biased dispersal and spatial heterogeneity affect landscape resistance to gene flow in fisher

Genetic connectivity results from the dispersal and reproduction of individuals across landscapes. Mammalian populations frequently exhibit sex-biased dispersal, but this factor has rarely been addressed in individual-based landscape genetics research. In this study, we evaluate the effects of sexbiased dispersal and landscape heterogeneity on genetic connectivity in a small and isolated population of fisher (Pekania pennanti). We genotyped 247 fisher samples collected across the southern Sierra Nevada Mountains of California. We tested for genetic evidence of sex-biased dispersal using sex-specific population structure and spatial autocorrelation analyses, and sex-biased dispersal tests of the assignment index, FST, and FIS. We developed resistance surfaces for eight landscape features hypothesized to affect gene flow and optimized each resistance surface independently by sex. Using multiple regression of distance matrices and an information-theoretic model selection approach, we fit models of genetic distance to landscape resistance distance separately by sex and geographic region. We found genetic evidence of sex-biased dispersal with significant differences in FST, FIS, and spatial autocorrelation between sexes. Optimal resistance values differed by sex, and model variables, fit, and parameter estimates varied substantially both between sexes and between geographic regions. We found a stronger relationship between landscape features and genetic distance for females, the philopatric sex, than the more widely dispersing males. Our results show that landscape features influencing gene flow differed by both sex and regional heterogeneity. Conducting analyses by sex and by region allowed for the identification of landscape genetics relationships not discernible when analyzed together. Our results show that failing to account for these factors can confound results and obscure relationships between landscape features and gene flow.


INTRODUCTION
Connectivity, defined as the ability of organisms to move within and among populations, is fundamental for long-term species persistence (Lowe and Allendorf 2010). For small populations vulnerable to stochastic events, demographic and genetic connectivity is essential for maintaining population viability (Gilpin and Soule 1986, Pierson et al. 2015, Benson et al. 2016). Genetic connectivity results from successful dispersal and reproduction of individuals across a landscape. Therefore, understanding the dispersal characteristics of a species and gene flow is inherently intertwined. Estimating dispersal through direct methods (e.g., mark-recapture, telemetry) is difficult for many species as it can be cost prohibitive over large landscapes, especially for species that are rare or difficult to detect. Direct methods also cannot easily distinguish between movement and gene flow, in which movement is followed by successful reproduction. Genetic methods for assessing dispersal (Slatkin 1987) coupled with non-invasive genetic sampling provide a powerful, cost-effective tool for sampling individuals over large areas and improving our understanding of dispersal and gene flow.
Landscape features can influence the ability of an individual to disperse through the landscape, and consequently, landscape features influencing dispersal will also influence gene flow. Previous studies have shown that gene flow can be influenced by a wide variety of factors including abiotic factors such as topography (Blair et al. 2013), climate (Schwartz et al. 2009), or anthropogenic features (Epps et al. 2013) and biotic factors such as land cover type (Cushman et al. 2006) or distribution of predators . Analysis of the influence of landscape features on gene flow can be either population-based, where individuals are grouped into demes and genetic differentiation is estimated among these demes, or individual-based, where analyses are conducted using individual genotypes without a priori definition of populations allowing for the assessment of fine-scale genetic structure across the landscape (Manel et al. 2003).
Sex-biased dispersal is a well-documented characteristic of many wildlife populations (Pusey 1987), and males and females may not respond to landscape features the same way during dispersal. In mammals, dispersal is often male-biased (Dobson 1982) although there are exceptions to this generalization (Favre et al. 1997). There are many hypotheses as to why sex-biased dispersal occurs, including resource competition (Greenwood 1980), inbreeding avoidance (Pusey 1987), and local-mate competition (Dobson 1982). Understanding sex-biased dispersal is important in conserving population connectivity, as factors important to connectivity for the dispersing sex may differ from those important to the philopatric sex and effective conservation strategies must address these differences.
The relative influence of particular landscape features on dispersal can also vary spatially. This may be due to biotic factors such as variation in population density, where individuals in dense populations may be more likely to disperse across landscape features than individuals from low-density populations (Matthysen 2005), or interspecific interactions where the presence of predators or competitors may influence dispersal behavior (Rundle and Nosil 2005). Alternatively, this spatial variation can be due to the differential availability of a feature on the landscape. Short Bull et al. (2011) identified the importance of replicating study areas in landscape genetics finding gene flow in American black bears differed among study areas due to variability of landscape features in each area. For example, in some of their populations, elevation was constant and had no relationship with gene flow, while in other areas elevation varied substantially and was found important for explaining gene flow. This conclusion has been found in a number of subsequent landscape genetic analyses (Moore et al. 2011, Trumbo et al. 2013. In this study, we examine how sex-biased dispersal and landscape heterogeneity affect gene flow in a small and isolated population of fisher (Pekania pennanti) in the southern Sierra Nevada Mountains of California. Identifying landscape elements structuring genetic connectivity and dispersal in this fisher population has important conservation implications. Due to its small size, estimated at <300 adults , and long-term genetic isolation (Knaus et al. 2011, Tucker et al. 2012, there are acute conservation concerns regarding the long-term viability of this population which has been considered for listing as a Threatened species under both State and Federal Endangered Species Acts (U.S. Fish and Wildlife Service 2014, California Department of Fish and Wildlife 2015).
The southern Sierra Nevada fisher population provides a useful study system to address questions of how sex-biased dispersal and landscape heterogeneity relate to genetic connectivity for a number of reasons. Sex-biased dispersal has been documented in other fisher populations (Powell 1993, Aubry et al. 2005, Matthews et al. 2013) and related mustelid species (Vangen et al. 2001, Zalewski et al. 2009) and so we hypothesized that this is likely also a feature of the Sierra Nevada fisher population. Secondly, this population occupies a large and diverse landscape with substantial spatial variation across the study area and a previous study has identified three genetic subpopulations that can be used to define groups for spatial replication of landscape genetic analyses (Tucker et al. 2014). Thirdly, fisher habitat has been well studied within California (Zielinski et al. 2004a, b, Davis et al. 2007) and elsewhere in their range genetic studies have identified landscape features that have been influential for fisher gene flow and provide a diverse suite of candidate variables from which to form landscape genetic hypotheses (Carr et al. 2007, Garroway et al. 2011, Hapeman et al. 2011).
Here, we use an individual-based landscape genetics approach to assess the importance of landscape features in influencing patterns of gene flow, and whether these landscape features vary by sex. Using samples from across the southern Sierra Nevada fisher population, we address the following questions: (1) Is there genetic evidence of sex-biased dispersal? (2) What landscape features influence genetic connectivity? (3) Does the influence of landscape features on gene flow vary by sex? and (4) Does the influence of landscape features on gene flow vary spatially across the study area? We investigate these questions by first testing for genetic evidence of sex-biased dispersal and then employ a resistance modeling approach in which we generate resistance surfaces representing our hypothesized relationships between landscape features and gene flow. This study contributes to the growing body of research that has detected sex-biased dispersal using a landscape genetics approach (Coulon et al. 2004, Shafer et al. 2012, Amos et al. 2014) but adds the unique elements of optimizing resistance values by sex, employing a sex-specific information-theoretic model selection approach, and comparing sex-specific to sex-combined analyses. This last step is important in that it allowed us to discern the potential problems that could arise in landscape genetics studies if you fail to account for sex bias when it is a characteristic of the species under study.

Genetic sampling
Fisher hair samples were collected from across the southern Sierra Nevada fisher population from 2006 to 2009 in conjunction with the U.S. Forest Service Sierra Nevada Carnivore Monitoring Program , 2017. The study area is~12,240 km 2 and encompasses the extant fisher population in the southern Sierra Nevada. Samples were collected from 223 survey units established across the elevational range of fisher (~1000-3400 m) on a systematic survey grid co-located with the Forest Inventory and Analysis sampling grid (Roesch and Reams 1999). Sample units consisted of arrays of six detection stations covering~0.8 km 2 with a mean minimum distance between sample units of 4.1 km. Each detection station was comprised of a baited track plate box with a barbedwire hair snare. An average of 140 units were sampled from May to October of each year. Additionally, in areas potentially occupied by fisher but not covered by the monitoring sampling frame (primarily Yosemite, Sequoia, and Kings Canyon National Parks), we opportunistically established hair snares along established roads or trails (at least 50 m from the road or trail) to fill in any potential spatial gaps in sample distribution. Details regarding study area and sampling methods are detailed in Zielinski et al. (2013) and Tucker et al. (2014). All hair samples from stations that detected fisher via tracks were genetically analyzed. We genotyped 247 fisher samples using 10 microsatellites (Dallas and Piertney 1998, Davis and Strobeck 1998, Duffy et al. 1998) and the y-linked marker DBY-3 (Hedmark et al. 2004). A y-locus control sample from a known male fisher was used to distinguish failed PCR from true negatives (females). Individuals captured multiple times were assigned the coordinates of their initial detection.

Sex-biased dispersal
We first assessed sex-biased dispersal by conducting a genetic structure analysis of males and females independently using the package GENE-LAND 4.0.3 (Guillot et al. 2005) in R 3.1.2. (R Development Core Team 2014) using the spatial uncorrelated allele frequencies model with the following parameters: K = 1-10; Marcov Chain Monte Carlo (MCMC) iterations = 500,000; thinning = 100; spatial uncertainty = 2000 m. Each model was replicated 20 times to evaluate consistency. To compare how strongly males and females assigned to a subpopulation, we calculated the mean and standard deviation of the probability of assignment of each individual to their most likely subpopulation. We would expect that with sex-biased dispersal, the philopatric sex would assign with greater certainty to the subpopulation of origin than the dispersing sex.
Next, we assigned individuals to subpopulations using the results of the sex-specific genetic structure analysis in order to test for differences between males and females in three different population-based genetic metrics (Goudet 2001): the assignment index (AI; Goudet et al. 2002), F ST (Weir and Cockerham 1984), and F IS . The AI estimates the probability of each individual's genotype originating within its geographic subpopulation of capture such that the dispersing sex should have a lower mean and higher variance of AI than the philopatric sex (Goudet et al. 2002). For the F ST test, we expect the philopatric sex to have significantly higher F ST values than the dispersing sex. F ST and AI tests were conducted across all subpopulation comparisons as well as for each pair of adjacent subpopulations. Significance of sex-biased dispersal tests for AI and F ST was determined using 10,000 randomizations using FSTAT 2.9.3 (Goudet 2001). We also estimated F IS by sex for each subpopulation with the expectation that with sex-biased dispersal, the dispersing sex should have a positive F IS (heterozygote deficiency) caused by a Wahlund effect because the dispersers' genotypes originated in a different subpopulation from where they were sampled (Wahlund 1928). F IS was calculated in GENEPOP 4.2 with significance determined using exact tests for heterozygote deficiency (Rousset 2008). Because our sampling was moderately male-biased, and there is debate as to whether this sampling bias may significantly influence results (Bekkevold et al. 2004, Hammond et al. 2006, we examined the effect of unequal sample sizes between males and females. We created 20 subsampled datasets in which males were randomly sampled without replacement to an equivalent sample size as females and repeated our analyses (Yannic et al. 2012). As other factors can produce a positive F IS (null alleles, positive assortative mating), we conducted a linear regression of F IS and F ST across all loci as a positive relationship between the two is indicative of a Wahlund effect (Waples 2014).
Lastly, we also employed an individual-based method to test for sex-biased dispersal using multi-locus spatial autocorrelation analyses in GenAlEx 6.5 Smouse 2006, Banks andPeakall 2012). We calculated the genetic autocorrelation coefficient (r) for each distance class and calculated 95% CIs for each estimate of r using 999 bootstraps from within each distance class. We defined the 95% CI around the null hypothesis (r = 0) by 1000 random permutations of the genetic distance matrix. The detection of spatial autocorrelation is influenced by distance class size (Peakall et al. 2003), and there is a tradeoff between spatial resolution and sample size per class. Therefore, we conducted analyses using 6-and 12-km increments to derive two different sets of distance classes. The 6-km size reflects the average female fisher dispersal distance in the southern Sierra Nevada (5.8 km; Sweitzer et al. 2015) and the optimal distance class size has been found to correspond with the mean dispersal distance of the philopatric sex (Banks and Peakall 2012). The 12-km size allows for increased sample size within classes but is still well below the maximum dispersal distance for female fisher in the Sierra Nevada (24.5 km, Sweitzer et al. 2015). We used two statistical criteria to assess the significance of the results from the spatial autocorrelation analysis: (1) r values for the philopatric sex were significantly greater than for the dispersing sex with 95% bootstrap CIs that do not overlap and (2) a heterogeneity test between sexes for each distance class using a squared paired-sample t-test (T2; Banks and Peakall 2012).

Landscape genetics-resistance modeling
We identified eight landscape features that we had a priori hypotheses regarding their effect on gene flow (Table 1; Appendix S1). We hypothesized that six landscape features were resistant to genetic connectivity including major water bodies (rivers/lakes; Wisely et al. 2004, Garroway et al. 2011, roads (Garroway et al. 2011), steep slopes (Jordan 2007), open areas (Powell 1993), moderate-or high-severity fire areas , and both very low and very high elevations (Davis et al. 2007. We also hypothesized that large tree size and high canopy cover would facilitate gene flow (Zielinski et al. 2004a, Purcell et al. 2009) and included two variables reflecting these hypotheses. We defined these variables using the California Wildlife Habitat Relationship system (Mayer and Laudenslayer 1988) and included potentially suitable forest types with large trees (>28 cm) and dense canopy cover (>60%). The two variables differ in the presence of the Sierra Mixed Conifer (SMC) habitat type, where its exclusion from the dense forest variable restricted the elevation range of forested pixels to the core of the occupied fisher range in the Sierra Nevada (1400-2300 m; Spencer et al. 2011), whereas its inclusion in the second variable (SMC forest) expanded forested pixels to the margins of the fisher elevation range. These higher-elevation areas may not be optimal habitat for fisher occupancy, but we hypothesized that they may play an important role in fisher gene flow.
Resistance layers were generated in ArcGIS 10 using a pixel size of 100 m. Variables were coded as resistance surfaces with features resistant to gene flow assigned a value greater than one. For example, for the two forest variables hypothesized to be conductive to gene flow (dense forest, SMC forest), all forested pixels were set = 1 (no resistance) and all non-forested pixels >1 (resistance). Conversely for water, hypothesized to be resistant to gene flow, all water pixels were assigned a value >1 (resistance) and non-water pixels = 1 (no resistance). As map edges have been shown to bias resistance values (Koen et al. 2010), we buffered locations by a minimum of 12 km. Resistance distance calculations require only one individual per pixel so if two individuals occupied the same pixel, we moved one of them~100 m in a random direction (McRae and Shah 2009). We assumed that this had a negligible effect because our sampling involved baited stations and scent lure which are thought to extend the sampling area hundreds of meters beyond the sampling device Mori 2001, Schlexer 2008).
We used circuit theory to estimate resistance distances between all pairs of individuals using CIRCUITSCAPE (McRae 2006, McRae et al. 2008). CIRCUITSCAPE incorporates multiple Notes: Each variable was optimized as both continuous and categorical variables over a range of maximum resistance values (2-100). The optimum resistance in the final column lists the data type (categorical or continuous) and resistance value that maximized the correlation between genetic distance and resistance distance. Species composition of CWHR types is detailed in Mayer and Laudenslayer (1988). alternative paths, accounts for effects of path size, and is directly related to random walks, which more realistically approximates dispersal behavior than the alternative least cost path approach, which assumes that dispersers have perfect knowledge of a single optimal path through the landscape (McRae et al. 2008, Pinto and Keitt 2009, Sawyer et al. 2011. CIRCUITS-CAPE analyses were conducted in pairwise mode, with individual locations set as focal nodes, and each pixel connected to its eight neighboring pixels. We established a null model representing Euclidean distance by running CIR-CUITSCAPE with a homogeneous resistance surface, with all pixels = 1. This null model has been established as the appropriate surrogate for geographic distance between individuals for circuit theory analyses instead of pure Euclidean distance (Koen et al. 2010, Amos et al. 2012. All references to Euclidean distance in our analyses refer to this null model. We calculated pairwise genetic distance between individuals using two metrics to compare results: the proportion of shared alleles (Dps; Bowcock et al. 1994), calculated using MSA 4.05 (Dieringer and Schlotterer 2003), and Rousset's a, calculated using SPAGeDi (Hardy and Vekemans 2002).

Optimizing resistance values
Assigning resistance values to landscape features is a well-acknowledged challenge in landscape genetics (Rayfield et al. 2010, Spear et al. 2010. To address this uncertainty, we evaluated each landscape feature over a wide range of resistance values from 2 to 100 (Table 1) and assessed the correlation strength with genetic distance using partial Mantel correlation coefficients (r pm ) (Mantel 1967, Smouse et al. 1986). Thematic resolution can also affect the nature and strength of landscape genetic relationships ) and so we also evaluated landscape features both as continuous and as categorical variables. The partial Mantel test has been found to be highly sensitive to changes in underlying landscape features (Kierepka and Latch 2015) but there is an ongoing debate regarding their use due to potential bias in statistical significance tests and inflated type 1 error rates (Legendre and Fortin 2010, Graves et al. 2013, Guillot et al. 2013, Legendre et al. 2015. However, our optimization approach relied solely on relative comparisons of the partial Mantel correlation coefficient, not the statistical test. This optimization method has been repeatedly found effective in identifying landscape features influencing gene flow, especially when combined with a comparative model selection approach (Wasserman et al. 2013, Castillo et al. 2014, Mateo-S anchez et al. 2015. We conducted resistance value optimization for males and females separately and with both sexes combined using the Ecodist package in R (Goslee and Urban 2007).
We defined optimum univariate resistance as the asymptote of the curve (rate of change <5%) of the partial Mantel r against resistance values, and employed two different optimization approaches: (1) independent optimization, where the optimum resistance value was defined individually for each sex, and (2) combined optimization, where the optimum resistance value was based on the overall strongest relationship with genetic distance for either sex. Previous research has identified the potential importance of also conducting multivariate resistance value optimization (Shirk et al. 2010). Therefore, we used these optimized univariate resistance values to create a candidate suite of multivariate models and conducted a model selection approach to identify an optimal multivariate resistance surface (details below).

Model selection and evaluation
We used optimized resistance values to fit linear models of landscape resistance to genetic distance with statistical significance assessed using multiple regression of distance matrices (MRDM). Multiple regression of distance matrices is the multiple regression of a response distance matrix (pairwise genetic distance) against two or more environmental variables (pairwise resistance distance) with significance testing performed by random permutations of the response distance matrix (Lichstein 2006, i.e., DIST GENETIC~D IST RESISTANCE1 + DIST RESISTANCE2 + Á Á Á + DIST RESISTANCEn ). Multiple regression of distance matrices has been found to have high power and low type 1 error for landscape genetic analyses (Balkenhol et al. 2009).
To test our first hypothesis that landscape features influence males and females differently, we fit an identical suite of models for each of three groups: all individuals, males only, and females only. To test our second hypothesis that landscape features important to gene flow vary spatially, we divided the data into two geographic regions consisting of adjacent subpopulation pairs based on the subpopulation structure identified in Tucker et al. (2014; North-Central and Central-South; Fig. 1). We chose to divide the data into two regions vs. three individual subpopulations as the latter would only provide insight into how the landscape influences gene flow within each subpopulation, whereas the regional division also speaks to what landscape features may be limiting gene flow between subpopulations. Subdividing in this manner also excludes the largest pairwise geographic distances (north-south), which on average (168 km) far exceeded the maximum reported dispersal distances for fisher both in the southern Sierra Nevada (36.2 km, Sweitzer et al. 2015) and in all of North America (~107 km, York 1996). Including such long-distance pairs may overestimate the importance of rare long-distance dispersal over more frequent short-distance dispersal (Parks et al. 2013). All individuals in the central subpopulation were included in both the north-central and central-south regions. However, due to the nature of the pairwise data used in this analysis, this did not create pseudo-replication as each region contained unique pairs of individuals between subpopulations. The combination of partitioning both by sex and by region resulted in nine separate analyses (all individuals, males, and females X all regions, north-central region, and central-south region). We then used Akaike's information criterion (AIC) model selection to identify the model best supported by genetic data.
Prior to creating a candidate model set, we assessed multicollinearity between variables by calculating Pearson's correlations and variance inflation factors (VIF) for a global linear model containing all variables. During this initial assessment, we found very high collinearity between many variables and between each landscape variable and Euclidean distance (Pearson's correlation >0.95, VIF >100). We discerned that this collinearity was generated in part because each landscape resistance distance matrix includes Euclidean distance. To address this problem, we subtracted pairwise Euclidean distances (null resistance surface) from landscape resistance distances prior to fitting MRDM models. Subtracting the null resistance surface from the landscape resistance surface enables analysis of how much more or less important a pixel is than expected under the null model (B. McRae, personal communication). We then re-calculated Pearson's correlations and variables with a correlation >0.80 were excluded from being in the same candidate model.
We created a set of 23 candidate models (four univariate and 19 multivariate) shaped by three overarching hypotheses that fisher gene flow (1) is facilitated by mid-elevation dense forest found to be strongly associated with fisher occupancy and denning habitat , (2) is facilitated by forested areas across a broad elevation range including high-elevation forests, and (3) is impeded by open areas such as large water bodies, roads, and burned areas (Appendix S2). Because the optimum resistance value varied, we standardized variables using z-transformations to facilitate comparison of parameter estimates and then fit MRDM models of genetic distance to landscape resistance using 10,000 permutations for significance tests. Models were again assessed for multicollinearity, and models with a variable VIF >5 were dropped from the candidate model set.
Models were ranked using second-order AIC c values, and models with uninformative parameters (additional parameter that did not improve the AIC score by at least 2) were excluded from further consideration (Burnham andAnderson 2002, Arnold 2010). We report parameter estimates from the top model-rather than modelaveraged coefficients due to concerns regarding model averaging and potential collinearity among variables (Cade 2015). We considered variables informative if the 95% confidence intervals (CIs) of their slope coefficients excluded zero. We report 95% CIs for parameter estimates, rather than the 85% CIs suggested for AIC-based approaches (Arnold 2010), to be more conservative in our assessment of potential landscape effects. There have been concerns regarding the use of AIC with MRDM, as pairwise data may violate the assumption-independent data points (Burnham and Anderson 2002). However, this approach has been regularly used in landscape genetics analyses in assessing model performance because the potential error is equivalent for each model, and consequently does not affect model ranking (Garroway et al. 2008, Richardson 2012, Engler et al. 2014). We conducted ❖ www.esajournals.org model fitting and selection in R using the following packages: GenABEL for z-transformations (Aulchenko et al. 2007), car for VIF calculation (Fox and Weisberg 2011), ecodist for MRDM (Goslee and Urban 2007), and MuMIn for model selection and averaging (Barton 2013). We investigated the effect of unequal sample size between sexes by generating 20 random male subsets, equivalent to female sample sizes, and repeating the model selection procedure.
We used the coefficients from the top model equations to generate a single optimum multivariate resistance surface for each sex and region. If our approach was effective, we would expect to see an increased r pm correlation for the multivariate surfaces compared to univariate surfaces. We used the raster calculator in the Spatial Analyst extension of ArcGIS 10 to multiply the untransformed parameter coefficient by its corresponding resistance surface and then summed all surfaces to generate a single multivariate resistance surface (Garroway et al. 2011, Row et al. 2015. Finally, using CIRCUITSCAPE, we generated pairwise resistance distances for these multivariate surfaces and calculated their r pm correlation with genetic distance.

RESULTS
Genotyping resulted in a dataset of 124 individuals comprised 69 males, 48 females, and seven individuals for which we were not able to determine sex. Sample sizes by subpopulation were as follows: north n = 42 (16 females, 22 males, four unknown), central n = 31 (12 females, 18 males, one unknown), and south n = 51 (20 females, 29 males, two unknown) (Fig. 1). We found limited movements in individuals recaptured multiple times with most recaptures occurring either within the same sample unit or at an adjacent sample unit. Genetic distance metrics were highly correlated (Dps~Rousset a: Mantel r = 0.94, P = 0.001), and so for clarity, we only report the results of the analyses using Dps.

Sex-biased dispersal
We found evidence of sex-biased dispersal in both population-and individual-based analyses. GENELAND identified three genetic clusters (K = 3) for both males and females with similar boundaries between subpopulations for each sex. Females assigned more strongly to subpopulations (mean prob assign female = 0.98 [standard deviation (SD) = 0.04], male = 0.87 [SD = 0.09], P < 0.0001, t-test) and had higher F ST values than males ( Table 2). The sex-biased dispersal F ST test was significant across all subpopulation comparisons and the north-central region, but non-significant in the central-south region (Table 3).
We also found an excess of homozygotes for males, but not females, as is expected with malebiased dispersal (F IS-north : male = 0.053, female = 0.006/F IS-central : male = 0.043, female = À0.005/F ISsouth : male = 0.043, female = 0.033) but these differences were not statistically significant. Linear regression of locus-specific plots of F IS and F ST showed a positive correlation in the northcentral region (slope = 0.58, r = 0.60, P = 0.07) as ❖ www.esajournals.org expected with a Wahlund effect, but not in the central-south region (slope = À0.063, r = 0.05, P = 0.89). All tests for the mean and variance AI were non-significant (Table 3). Subsampling to balance the sample size of males and females did not produce a meaningful difference in results. Spatial autocorrelation analyses for both regions at 6 and 12 km distance increments found significant autocorrelation in females in the 24 km distance class and r female > r male (6 km: Appendix S3, 12 km: Fig. 2C-F). Results of sex-biased dispersal tests varied depending on the increment size. For the 12-km increment size, tests were positive for sex-biased dispersal in both regions at the 24 km distance class with no overlap in the 95% bootstrap CIs and significant T2 statistics ( Fig. 2A, B). For the 6-km increment size, sexbiased dispersal tests were also significant in the central-south region; however, while the northcentral had r female > r male , there was overlap in the bootstrap CIs and a non-significant T2 statistic according to the recommended P < 0.01 significance cutoff for the T2 heterogeneity tests (Banks and Peakall 2012, Appendix S3).

Landscape genetics-resistance modeling
Optimizing resistance values.-Results of landscape resistance analyses support that the influence of landscape features on gene flow varies by sex. When optimizing resistance values for each variable, we found that the partial correlations (r pm ) of genetic distance and resistance distance were markedly different for males and females for the majority of variables (Fig. 3). When analyzing the sexes together, r pm was generally intermediate between the values for each sex. Optimum resistance values differed considerably between variables ranging from 2 (slope) to 100 (dense forest; Table 1), but were generally similar in both the independent and combined optimization approaches. We did not find any significant difference in the model selection results between these two approaches, and for simplicity, we only present the results of the combined optimization approach.
Model selection and evaluation.-Model selection also varied substantially between sexes in terms of variables in the top models, magnitude of parameter estimates (b), and model fit (R 2 , r pm ). Full AIC model selection results are provided in Appendix S4. Across all subpopulations, the top female model indicated that gene flow was impeded by large bodies of water and roads, but facilitated by dense forests. For males, the top model was univariate, comprising solely the SMC forest variable, indicating that for males gene flow was facilitated by SMC forest (Table 4). There were no variables in common between top models for males and females. With sexes combined, the top model was a mix of variables in the top models for each sex, but with reduced parameter coefficients and model fit compared to the sex-specific analyses. Explanatory power of models for males, and for the sexes combined, was low (R 2 = 0.05-0.10).
Results also supported spatial variation in the influence of landscape features of gene flow as top models varied by region (Table 4). Respective forest variables for both sexes (females = dense forest, males = SMC forest) were consistent across regions, but all other variables in the top models differed. Gene flow in females was most strongly influenced by water in the north-central region and by roads in the south-central region (Table 4). For males, the north-central top model found gene flow to be influenced by SMC forest and roads, but this model had high uncertainty and poor fit (w i = 0.51, R 2 = 0.05). In the central-south region, SMC forest was the only variable in the top model for males. Taking subsets of male data to an equivalent sample size as females produced generally consistent results as the full male dataset although a minority of the male subsets included the alternate variables of either elevation or openings instead of SMC forest in the top model. This is likely due to a high correlation between these variables (>0.70) as they are all comprised of broadly distributed areas of low resistance across a wide elevational range. In the north-central region, the male subsets had poor fit and high model uncertainty resulting in variation in the covariates in the top model, which is indicative of the low explanatory power of landscape features on male gene flow in this region. Notably in both regions, some male subsets found Euclidean distance as the top model (Appendix S5). Females had better model fit than either males alone or the sexes combined, and the composite multivariate resistance surfaces for females had increased partial Mantel correlations compared to univariate surfaces. The female north-central Correlograms for sex-and region-specific spatial autocorrelation analyses. X-axis: distance class endpoint in 12-km increments; y-axis: genetic correlation coefficient (r) for females (solid lines) and males (dashed lines). (A, B) comparison between sexes by region with 95% bootstrap confidence intervals with asterisks Ã indicating statistically significant sex-biased dispersal tests. (C-F) sex-and region-specific results with permutation tests (dotted lines) illustrating the null model of no spatial autocorrelation: C = females north-central; D = males north-central; E = females central-south; and F = males central-south. regional model had the best fit and highest correlation for the composite multivariate resistance surface of any sex or region (R 2 = 0.25, r pm = 0.31). In the north-central region, the model for sexes combined had equivalent variables as the female-only model. However, across all subpopulations and in the central-south region, the models for the sexes combined were a mix of variables from the male and female models. For all analyses, parameter coefficients with the sexes combined were reduced compared to the same variable in the sex-specific analysis. Top models for all sexes and regions were statistically significant at P < 0.05.

Maximum Resistance
ParƟal Mantel CorrelaƟon (r pm ) Fig. 3. Plots of the partial Mantel correlation value (r pm ) of pairwise genetic distance~resistance distance, controlling for the effects of Euclidean distance, for a range of maximum resistance values for each variable (2-100). Correlations are plotted for all individuals (purple circles), females (red triangles), and males (blue squares) independently.

DISCUSSION
Our landscape resistance modeling results show that both sex-biased dispersal and landscape heterogeneity can influence the association between landscape features and gene flow. Conducting resistance modeling separately by sex and region allowed for the identification of landscape genetics relationships not discernible when analyzed together, and failing to account for these factors may confound results and obscure relationships between landscape features and gene flow.
For analyses conducted with the sexes combined, we found the landscape relationships for one sex often masked one or more variables important to gene flow for the other sex. This was due to that sex either having a larger sample size or having a stronger relationship to genetic distance. With the sexes combined, parameter estimates of variables were lower than those in the sex-specific models, due to the averaging of opposing responses between the sexes. This issue was also reported by Amos et al. (2014) who detected strong sex-specific responses to habitat fragmentation when conducting analyses by sex but had inconclusive results when analyzing both sexes together.
We also found regional variation in landscape features influencing gene flow. Only the forest cover variables for each sex (females: dense forest, males: SMC forest) were consistent in the top models, with all other variables differing between regions. Regional differences were most pronounced when analyzing the sexes combined where there were no variables in common between the two regional models. Our results concur with other studies that have found spatial heterogeneity in landscape genetic relationships, and we echo the warning that it may not be appropriate to apply landscape genetic findings from a specific region to other populations within a species (Short Bull et al. 2011, Trumbo et al. 2013).

Sex-biased dispersal
Overall, the results of sex-biased dispersal tests provide support that males disperse farther than females. F ST and F IS values were consistent with male-biased dispersal with females having significantly greater F ST and stronger subpopulation assignment, F IS > 0 for males but not females, and significant male bias detected in spatial autocorrelation analyses. Notably while the F ST test was significant for the entire population, regionally it was only significant for the north-central pair. This is likely being driven by generally higher gene flow between the central and south groups as evidenced by comparatively low F ST values for both males and females in this region. There was no significant difference between sexes in either the mean or variance of AI. However, Goudet et al. (2002) found that except when the dispersal rate is very low (<10%), the F ST statistic is the most powerful and robust to changes in sampling and magnitude of sex bias for detecting sex-biased dispersal. Dispersal estimates for fisher primarily focus on juvenile dispersal and vary widely from a mean of 16.7 km-female/41.3 km-male in British Columbia (Weir and Corbould 2010) to 1.3 km-female/ 4.0 km-male in northern California (Matthews et al. 2013), and 5.8 km-female/9.8 km-male in the southern Sierra Nevada (Sweitzer et al. 2015). While less is known about adult dispersal, males have been found to move broadly during the spring breeding season (Sweitzer et al. 2015). Therefore, dispersal characteristics of fisher may result in low power using AI-based test statistics. Additionally, such mixed results of statistical testing have been reported in studies of species that otherwise have strong evidence for sexbiased dispersal due to difficulties in using biparentally inherited markers described in more detail below (Helfer et al. 2012).
There are inherent difficulties in detecting sexbiased dispersal using bi-parentally inherited markers and the timing of sampling. With nonoverlapping generations, microsatellite markers can detect sex-biased dispersal when individuals are sampled after dispersal, but before reproduction, when alleles from dispersers will be transmitted to its offspring. The presence of overlapping generations, more realistic of wild populations, extends the genetic signal of sex-biased dispersal as mothers and daughters cluster geographically, while male offspring disperse farther from their natal area (Blyton et al. 2015). Juvenile dispersal of fisher from their natal areas usually begins in mid-late winter (~10 months post-birth), with the majority of juvenile dispersal occurring from February to September (Lofroth et al. 2010, Sweitzer et al. 2015. Because our sampling period (June-October) overlaps this dispersal period and also spans multiple years, we could potentially sample the same individuals both pre-and postdispersal. However, the majority of recaptured individuals were detected in close spatial proximity to their initial detection and we did not recapture the same individuals in multiple subpopulations, indicating that the individuals sampled were not likely dispersing long distances over the course of the study.

Landscape genetics-resistance models
Our results indicate that the influence of landscape features on gene flow varies by sex and that landscape resistance analyses conducted for males and females independently are more biologically relevant for fisher in the southern Sierra Nevada than the analyses conducted with the sexes combined. We reached this conclusion based on (1) differences between sexes in correlation values during optimization of resistances surfaces, (2) differences between sexes in variables and parameter estimates of top models, and (3) improvement in model fit when analyzing the sexes independently.
In all regions, dense forest had a positive influence on genetic connectivity of females. This relationship is consistent with home range characteristics for fishers in the Sierra Nevada (Zielinski et al. 2004a, Purcell et al. 2009) and throughout their range (Aubry et al. 2013). This relationship is also congruent with fisher denning habitat models that have found female den sites primarily in dense mid-elevation mixed conifer and hardwood forests at a comparatively narrow elevational band compared with more broadly distributed foraging habitat (Spencer et al. 2015). We found an entirely different suite of variables in the top models of males. The SMC forest variable was present across all male top models instead of the dense forest variable found important in females, and model fit for males was much lower than for females in all analyses. The comparatively poor model fit and low partial Mantel correlations for males imply that landscape features create relatively little resistance to gene flow and that the landscape overall has much less influence on male dispersal than on female dispersal.
Our results also indicate that the influence of landscape features on gene flow varies spatially. While the forest variables were consistent for each sex across regions, all other variables differed between regions. In the north-central region, water had a strong negative influence on genetic connectivity of females, but not males. The negative influence of water on female gene flow is likely attributable to the Kings River Canyon, which has long been hypothesized to limit fisher dispersal (Wisely et al. 2004, Jordan 2007, Tucker et al. 2014) and may also reflect fine-scale genetic structuring among several small clusters separated by the San Joaquin River within the northern subpopulation (Tucker et al. 2014). Notably, the water variable was absent from all male top models. Considering the magnitude of the Kings River Canyon, it is surprising that it does not appear to create any impediment to male dispersal. In the central-south region, roads rather than water had a negative influence on genetic connectivity for females. This relationship was unexpected as the road density and intensity of road use seem comparable across regions and even the most heavily used roads within the central-south region are primarily two lane mountain roads running through otherwise contiguous forested habitat. For males in the north-central region, roads were found to negatively influence gene flow, but not in the centralsouth region where the top model con.sisted solely of the SMC forest variable. Both regional models for males had very poor fit.
The difference between forest variables selected by males (SMC forest) and females (dense forest) is particularly intriguing. The SMC forest variable encompasses a higher elevational range than the dense forest variable (maximum elevation SMC forest = 2602 m, dense forest = 2223 m) and includes a much larger geographic area covering an additional 171,849 pixels (~1718 km 2 ) over the dense forest layer. We speculate that this indicates female gene flow is characterized by dispersal among high-quality habitat in the core elevation range for fisher occupancy in the Sierra (~1400-2300 m elevation; Spencer et al. 2011), while males disperse more widely and therefore associate with a more widespread land cover type.
Our findings indicate that female dispersal may be one of the factors limiting population expansion north of its current extent. The northern boundary of the current population aligns with two features shown to impede female fisher flow: the Merced River and the heavily traveled roads associated with Yosemite National Park. Habitat models indicate that there is suitable but unoccupied habitat north of the Merced River (Davis et al. 2007). There have been incidental sightings of fishers north of this river but no evidence of a breeding population (Chow 2009), and it is likely that these sightings were dispersing males that failed to find mates. Consequently, northward population expansion, a goal identified in a recent fisher Conservation Strategy (Spencer et al. 2016), may require assisted translocation of females across these potential dispersal barriers.
Previous research has shown a temporal delay between genetic data and landscape features, making it difficult to know whether to attribute current genetic patterns to contemporary or historical landscape characteristics. However, this lag time has been found to be relatively short when using Dps and Mantel correlations in comparison with population-based approaches using F ST . The temporal delay is also smaller for species with longer dispersal distances (>10 km, Landguth et al. 2010). Considering the relatively long dispersal ability for fisher reported for the southern Sierra Nevada (36.2 km, Sweitzer et al. 2015) and elsewhere in their range (~100 km, York 1996), our results are likely attributable to recent (within the last few decades) rather than historical landscape conditions.
There is an important distinction between population-based analyses (F ST , F IS , AI, population assignment) and individual-based analyses (spatial autocorrelation and landscape resistance modeling). This distinction is critical when considering the effect of dispersers in each approach. In a population-based approach, one individual dispersing between subpopulations may only have a small effect on population-based metrics like F ST , but in an individual-based approach, one disperser will affect pairwise values to all other individuals both within and between subpopulations. Consequently, population-based analyses will only reflect the influence of landscape features between subpopulations, whereas individual-based analyses reflect the influence of landscape features both between and within subpopulations (Broquet et al. 2006, Segelbacher et al. 2010. We acknowledge that the sample sizes used in this study are small, which decreases power of statistical analyses. In landscape genetics analyses, reducing the number of individuals does not affect the accuracy of the estimate of the partial Mantel correlation, but does increase the uncertainty in the estimate (Landguth et al. 2012).
However, the issue of sample size has never been directly addressed using MRDM and model selection. The sample sizes used in this study are similar to other individual-based landscape genetic analyses detecting significant relationships between genetic distance and landscape features (Cushman et al. 2006, Shafer et al. 2012, Wasserman et al. 2013. Considering this population has been estimated to contain <300 adults , the 124 individuals analyzed in this study represent a substantial proportion (≤41%) of the extant population.

CONCLUSIONS
Careful consideration of the potential for sexbiased dispersal and landscape heterogeneity should be undertaken prior to conducting landscape genetic analyses as these factors can strongly influence results. Our findings have considerable conservation implications in showing that failing to account for sex-biased dispersal in landscape genetics analyses may result in omitting or misidentifying landscape features important for genetic connectivity.
This study indicates that management actions with the goal of conserving or enhancing population connectivity need to consider that males and females may each have a different suite of habitat features important for connectivity. For species with male-biased dispersal, population expansion is likely mediated by female's ability to disperse to new habitat areas. Consequently, identifying and conserving landscape features important for female dispersal may be more important than for more widely dispersing males.