Comparison and parallel implementation of alternative moving-window metrics of the connectivity of protected areas across large landscapes

A variety of metrics can be used to measure connectivity of protected areas. Assumptions about animal movement and mortality vary among metrics. There is a need to better understand what to use and why, and how much conclusions depend on the choice of metric. We compare selected raster-based moving-window metrics for assessing the connectivity of protected areas to natural habitat in the surrounding area, and develop tools to facilitate calculation of these metrics for large landscapes. We developed parallel implementations of distance-weighted sum and Spatial Absorbing Markov Chain methods in R packages to improve their useability for large landscapes. We investigated correlations among metrics for Canadian protected areas, varying background mortality, cost of movement, mean displacement, dispersal kernel shape, distance measure used, and the treatment of natural barriers such as water, ice, and steep slopes. At smaller spatial scales (2–5 km mean displacement), correlations among metric variants are high, suggesting that any of the metrics we investigated will give similar results and simple metrics will suffice. Differences among metrics are most evident at larger spatial scales (20–40 km mean displacement) in moderately disturbed regions. Assumptions about the impact of natural barriers have a large impact on outcomes. In some circumstances different metrics give similar results, and simple distance-weighted metrics likely suffice. At large spatial scales in moderately disturbed regions there is less agreement among metrics, implying that more detailed information about disperser distribution, behaviour, and mortality risk is required for assessing connectivity.


Introduction
Protected areas are considered by many to be the foundation supporting efforts to conserve biodiversity (Coetzee et al. 2014;Watson et al. 2014;Gray et al. 2016). However, for protected areas to be effective, they need to be connected in ecological networks that allow for ecological and evolutionary processes, such as gene flow, colonization, migration, and species range shifts (Saura et al. 2018). Governments around the world formally acknowledged the importance of connectivity through the adoption of Aichi Target 11 at the 2010 Conference of Parties for the convention on biological diversity (CBD), which calls for Parties to establish and expand "well connected systems of protected areas and other effective area-based conservation measures" (CBD 2010). A similar but more ambitious target was adopted at the 2022 Conference of the Parties to the CBD as part of a new Global Biodiversity Framework (CBD 2020).
The public and governments require practical and meaningful metrics for measuring progress towards achieving the elements of Aichi Target 11, including metrics to assess how connected protected areas are (Geldmann et al. 2021). Indicators can also be used to help identify priority candidate areas for inclusion in protected area networks to increase the connectivity and effectiveness of the network. To support those objectives, three indicators were considered by the Subsidiary Body on Scientific, Technical and Technological Advice of the CBD for measuring progress towards the connectivity component of Target 11 (UNEP-WCMC 2020): protected connected (ProtConn; Saura et al. 2017), protected area representativeness and connectedness (PARC; Drielsma et al. 2007; Biodiversity Indicators Partnership 2019), and Beyer's ecosystem intactness index (BEII; Beyer et al. 2020). ProtConn has several variants but in general measures the percentage of an ecoregion covered by connected protected areas. ProtConn considers only the size and proximity of protected areas to one another and does not factor in the condition of the intervening landscape (Naidoo and Brennan 2019). In contrast, PARC considers the condition of the intervening landscape in two ways: by determining whether protected areas are connected to protected and non-protected natural areas, and by measuring distances using a least-cost path based on the cost of land-cover features in the intervening landscape. Similarly, BEII measures intactness based on a combination of habitat loss, fragmentation, and degradation resulting from anthropogenic disturbance. BEII was designed to measure the intactness of ecoregions but can be adapted to focus on landscapes in and around protected areas.
These metrics and the others we consider require one or more input maps (i.e. resistance to movement, habitat quality and/or mortality) that can be defined in a variety of ways. Typically, species-agnostic approaches are "based on the assumption that natural terrestrial areas facilitate connectivity, and the higher the degree of human modification, the more the 'ecological flow' is restricted" (Marrec et al. 2020). To focus attention on similarities and differences among the metrics (rather than the inputs) in this study we extend that assumption to derive species-agnostic cost, habitat quality and mortality layers from Canadian human footprint and natural barrier information compiled by Pither et al. (2021).
Both PARC and BEII are raster-based metrics that can be considered part of a larger family of movingwindow distance-weighted sum metrics that measure average natural habitat quality in the vicinity of each pixel, giving higher weights to habitats that are nearer to the focal pixel and easier to access. This family of metrics also includes simpler moving window metrics of local habitat density (Miguet et al. 2017), and resistant kernel approaches (Compton et al. 2007;Cushman and Landguth 2012;Diniz et al. 2020). To understand relationships among these metrics, it is helpful to recognize that each is associated with a set of assumptions about animal movement and mortality (Table 1). We refer to this set of assumptions as the "model". Distance-weighted sum metrics are mathematically related to simple dispersal models (Hughes et al. 2015) (Table 1), and a criticism of these and other connectivity models is that they do not effectively distinguish or account for effects of landscape structure on movement behaviour and mortality (Diniz et al. 2020;Day et al. 2020;Fahrig et al. 2021). To help address this deficit, Fletcher et al. (2019) have suggested metrics derived from a spatial absorbing Markov chain (SAMC) model that explicitly accounts for dispersal cost and mortality. Several classes of metrics can be derived from the SAMC model; unconditional short-term metrics include the disperser distribution at time t and cumulative mortality distribution at time t, conditional short-term metrics track the fates of dispersers according to their starting locations, and long-term metrics describe outcomes after an infinite number of timesteps (Table 1 of Fletcher et al. 2019).
A challenge for distance-weighted sum and SAMC metrics is that computing them requires substantial processing time and, in some cases, memory. PARC and resistant kernels require calculation of leastcost paths which are particularly demanding in this regard. BEII and other Euclidean distance-weighted sum metrics involve simpler moving window calculations that are less demanding, and relatively easy to implement in parallel. Marx et al. (2020) calculate SAMC metrics by constructing a full transition probability matrix, which allows calculation of all classes of SAMC metrics, but requires a large amount of memory which limits utility for analysis of large landscapes.
All our selected metrics include habitat amount because there is strong evidence that increasing habitat amount is good for biodiversity (e.g. Watling et al. 2020). These metrics are conceptually related to one another (Table 1) and are also likely to be highly correlated with one another. Correlations among metrics that include both habitat configuration and amount are problematic if the goal is to distinguish effects of landscape configuration from effects of habitat amount (Fahrig 2017). However, our interest here is in distinguishing protected areas that are well connected to habitats of various qualities in surrounding areas from those that are not. Finding agreement among alternative models in this context can clarify when more information is not needed and simplify decision making.
In this paper, we investigate relationships among SAMC short-term disperser distribution metrics and comparable distance-weighted sum metrics. We also present a parallel iterative moving window approach for calculating SAMC unconditional short-term metrics (Fletcher et al. 2019) that is less memory intensive, faster, and useable for larger landscapes, and a parallel implementation of Euclidean distanceweighted sum metrics. To investigate correlations among metrics, we assessed the connectivity of Table 1 Mortality and movement assumptions associated with moving-window distance-weighted sum (DWS) and SAMC connectivity metrics Note that a variety of different metrics can be derived from the same set of assumptions about movement and mortality (e.g. Table 1 of Fletcher et al. 2019 natural areas within and around 5149 Canadian protected areas using a variety of distance-weighted sum (including PARC and BEII) and SAMC metric variants at a range of spatial scales. Canada provides a good real-world case study for assessing the behaviour of connectivity metrics due to the diversity of natural barriers such as glaciers, steep slopes, and large water bodies, and human impacts that span a full range from highly developed cities and farms, through moderately developed shared lands, to large wild areas (Locke et al. 2019;Hirsh-Pearson et al. 2022). This work helps clarify the similarities and differences among connectivity modelling methods that can be used to inform protected area planning and tracking of connectivity over space and time, while also addressing a more general need to explore the behaviour among connectivity modeling methods (Diniz et al. 2020).

Data sources
The starting point for our analysis is a recently developed, 300 m resolution, pan-Canadian movement cost layer, parameterized for terrestrial animals that are likely to use and move more successfully through natural land cover types (Pither et al. 2021), using more complete and accurate Canadian human footprint information (Poley et al. 2022;Hirsh-Pearson et al. 2022). This pan-Canadian layer includes both anthropogenic land cover features (e.g., roads, cities, croplands) and natural features considered to negatively affect the successful movement of animals (high elevations, steep slopes, glaciers, large lakes and rivers). Pither et al. (2021) used expert opinion to assign values of 0, 10, 100 and 1000 to 7 natural and 20 anthropogenic features, representing the cost of movement, and then assigned the highest cost from among the 27 land features to each raster pixel i. Four categories were selected because more would imply a level of precision higher than can be obtained from expert opinion. To derive a second cost layer that does not include the human footprint, we set the cost of anthropogenic features to zero, and assigned the highest cost from among the remaining 7 natural land features to each pixel i. We refer to the cost layers with and without the human footprint as C F i and C N i , respectively. The full list of land cover features and their associated costs are in supplementary Table S3, and additional information on the source of each land cover feature is in Pither et al. (2021). In this national analysis, we have not included effects of fencing and wildlife passages along highways because we lack a consistent national map of fence and passage locations.
No cost layer will work for all species; however, Pither et al. 2021 compared the current density map produced using their cost layer to GPS-collar data for moose (Alces alces), grey wolf (Canis lupus), Rocky Mountain elk (Cervus canadensis), and mountain caribou (Rangifer tarandus), and roadkill data for herpetofauna and moose, and found that individual mountain caribou, grey wolves, moose, and Rocky Mountain elk who traveled across larger areas were more associated with areas of high current density. Moose roadkills were more likely in areas of high current density, but herpetofauna roadkills were not. Together, these results suggest this particular cost layer is most relevant for large animals moving long distances.
We converted the two cost layers to habitat quality layers, one that includes the human disturbance footprint ( H F i ) and one that does not ( H N i ) , by rescaling cost layers with and without anthropogenic features ( C F i and C N i , respectively), and assuming that habitat quality and cost are inversely related to one another to obtain habitat quality maps that vary between 0 and 1: H i = 1 − C i ∕1000 . We also assumed that low quality, high cost areas produce fewer dispersers and are associated with higher mortality risk (Table 2). We recognize that movement cost, habitat quality, mortality, and production of dispersers are not generally perfectly correlated with one another (Diniz et al. 2020;Drake et al. 2021), and note that the metrics we consider can easily accommodate deviations from these assumptions where better information on movement cost, mortality and disperser production is available.
Protected and conserved areas were obtained from the Canadian Protected and Conserved Areas Database (Environment and Climate Change Canada 2020), which includes protected areas that meet IUCN categories I to VI and other effective area-based conservation measures (as defined by Pathway to Canada Target 1 2019), and from the United States protected areas database 2.0 (U.S. Geological Survey 2020) which includes areas defined as being "dedicated to the preservation of biological diversity and to other natural (including extraction), recreation and cultural uses, managed for these purposes through legal or other effective means" (U.S. Geological Survey 2022). To help interpret spatial variation in results we used a 2014 ecozone map of Canada (Agriculture and Agri-Food Canada 2014), which includes 15 terrestrial ecozones (Fig. 1). The degree of human modification varies substantially across Canada and among ecozones (Hirsh-Pearson et al. 2022).
Moving window distance-weighted sum metrics All distance-weighted sum metrics define the value at pixel i as a weighted sum of habitat quality H j , where each pixel j is within a specified distance of i, and this moving window calculation is repeated around each pixel on the landscape. These metrics all have a spatial scale or dispersal capacity parameter that determines the extent of the moving window, making them somewhat species-specific. They vary in the choice of distance metric (e.g. Euclidean or least cost path), the shape of the distance weighting function or kernel f (d ij ) (e.g. uniform or exponential), whether or not the distance weighting function is normalized (i.e. ∑ j f � d ij � = 1 ) and whether or not the focal habitat quality is also included in the distance-weighted sum. For example, PARC connectivity is defined as where d ij are least-cost path distances from pixel i to j (Drielsma et al. 2007), is a spatial scale parameter, and w is the raster cell resolution. Movement cost for each pixel is inversely proportional to habitat quality: We selected parameter values for high contrast (i.e. large differences between c max , and c min ) and low contrast (i.e. smaller differences between c max , and c min ) metric variants (Table 2) to examine sensitivity to the degree of contrast, and highlight similarities in outcomes from PARC and a Spatial Absorbing Markov Chain metric that includes both mortality and resistance to movement (Fig. 2). PARC is a least cost path exponential distance-weighted sum metric without a normalized kernel or focal habitat quality included in the sum. We note that PARC is also a form of resistant kernel (Compton et al. 2007;Cushman and Landguth 2012;Diniz et al. 2020), and that using a non-normalized distance-weighting function means that signal is lost where habitat quality is low. This signal loss is somewhat analogous to mortality but difficult to interpret biologically (Table 1).
In contrast, the Ecosystem Intactness Index (BEII; Beyer et al. 2020) is defined as where d ij are Euclidean distances, so I is a Euclidean exponential distance-weighted sum metric with a standardized kernel and focal habitat quality included in the sum. Note that multiplying by focal habitat quality H i is equivalent to assuming that mortality occurs after dispersal (Table 1). Following Beyer et al. we assume the scaling exponent z = 0.5.
For comparison to these distance-weighted sum metrics that have been considered by Subsidiary Body on Scientific, Technical and Technological Advice of the CBD, we also examined two simpler metric variants (Miguet et al. 2017). A simple standardized Euclidean exponential metric (E) is and a simple standardized Euclidean uniform metric where r is the radius of a circle. These simple metrics include no mortality or signal loss, and imply that dispersal probability depends only on Euclidean distance (Table 1). Uniform kernels are not biologically realistic, but they are often used in "scale of effect" studies (Miguet et al. 2017) so we include them here.
To tease apart effects of spatial scale from kernel shape (or type) we parameterize the uniform and exponential kernels in terms of mean displacement (Online Resource Fig S1). If f we assume that habitat quality H i is proportional to the density of organisms dispersing from location i then mean displacement is the average distance moved by dispersing organisms. For a circular uniform kernel, there  Table 2) assuming mean displacement of 20 km for the area around Mont-Tremblant Park (Fig. 1). PARC is more similar to SAMC B than SAMC R because the non-normalized distance-weighting function causes signal loss that is somewhat analogous to mortality (Table 1) is a simple relationship between mean displacement and radius: = 2r∕3 . For the exponential metrics (Eqs. 1-3) we define the spatial scale parameter in terms of mean displacement : = 2∕ . To limit computational cost, we set small exponential kernel values ( e − d ij <10 -6 ) to 0. For all metrics, we consider mean displacement values of 2, 5, 10, 20 and 40 km to reflect a range of natal and gap-crossing dispersal distances found across a range of mammals and birds (e.g. Albert et al. 2017). For reference, Online Resource Table S1 contains definitions of all parameters and variables defined in this paper, and all acronym definitions are included in Online Resource Table S2.

Spatial absorbing Markov chain metrics
As with the distance-weighted sum metrics, we begin by assuming that the initial distribution of organisms across the landscape is proportional to our index of habitat quality. Distance-weighted sum metrics can be interpreted as the distribution of organisms after a dispersal (Table 1), with a finite dispersal period implicitly determined by the mean displacement parameter. We construct an analogous SAMC (Fletcher et al. 2019) metric from the disperser distribution at time t, with mean displacement implicitly determined by the dispersal period. The SAMC model allows for spatial variation in the mortality of dispersers at each timestep, while assuming that the movement of dispersers from each pixel at each timestep depends on relative cost of dispersal among adjacent pixels (Table 1). Expressing the number of dispersers at location i and time t (B i,t ) as an iterative focal sum helps clarify the relationship to distanceweighted sum metrics: where m i is mortality at location i, c j is cost of movement at location j, and w is raster cell width (assuming 8-neighbour rule).
Mean displacement increases with the number of timesteps. We calibrate the time parameter t by setting movement cost to one and mortality to a low background level of m j = 0.00001, then measuring the mean displacement of organisms dispersing from a single focal pixel. We select times (40, 250, 1000, 4001, and 16,007) that correspond to mean displacements of 2, 5, 10, 20 and 40 km respectively on a 300 m resolution landscape. Without spatial variation in mortality or movement cost, SAMC distribution kernels have a shape that is intermediate between the exponential and uniform kernels (Online Resource Fig S1).
We investigate SAMC mortality + cost (B), mortality (M) and movement cost (R) metric variants. To examine effects of varying mortality (SAMC B and M) we set mortality inversely proportional to habitat quality: where quality varies between 0 and 1. Similarly, to examine effects of varying movement cost (SAMC B and R) we set cost inversely proportional to habitat quality (Eq. 2). Parameter values for high contrast (i.e. large differences between m max and m min or c max , and c min ) and low contrast (i.e. smaller differences between m max and m min or c max , and c min ) SAMC metric variants are given in Table 2. The high contrast scenarios yield displacement patterns that are more influenced by landscape pattern (Online Resource Fig S2). We note that it is also possible to construct SAMC metrics analogous to exponential Euclidean distance-weighted sum metrics by considering background mortality as settling probability, and projecting the long-term mortality distribution (Robert Fletcher pers. comm.), but that approach would not allow us to distinguish survivors from the dead on landscapes where mortality is heterogeneous.
Parallel implementation of SAMC and euclidean distance-weighted sum metrics An advantage of formulating SAMC unconditional short-term metrics ( we present an efficient parallel method for calculating SAMC unconditional short-term metrics that is less memory intensive, faster, and useable for larger landscapes than the samc R package methods (Marx et al. 2020). In our pfocal R package (https:// lands citech. github. io/ pfocal/), we also offer an efficient parallel analog to the focal function in the raster package (Hijmans 2022a) that can be used to speed calculation of Euclidean distance-weighted sum metrics. In our tests, computation time decreased 10 to 100-fold, depending on the size of the kernel, the number of available processors, and the size of the landscape (Figs S3 and S5). Our parallel SAMC implemention also requires substantially less memory, and can run on landscapes that are too large for the samc R package (Online Resource Fig S4). Memory requirements of our pfocal method are comparable to the focal function in the raster package (Hijmans 2022a), but the focal function in the terra package (Hijmans 2022b) requires less memory (Online Resource Fig  S6). See https:// github. com/ LandS ciTech/ LSTD-Conne ctivi ty-Paper/ blob/ main/ scrip ts/x_ examp leMet ricCa lcula tions.R for examples of metric calculations using the LSTDconnect and pfocal R packages.
A disadvantage of the iterative moving-window approach is that SAMC conditional metrics cannot be calculated because those metrics require keeping track of where dispersers come from. Long-term unconditional metrics could be approximated by iterating until the number of live dispersers is negligible, but we do not know if there is an advantage to that approach. The LSTDConnect package is intended to demonstrate the value of a parallel iterative movingwindow algorithm and support the analysis in this paper. Rather than continue to develop a fully functional alternative to the samc R package (Marx et al. 2020) we are investigating options for integrating iterative parallel moving-window methods into that package.

Omission of SAMC R and PARC metrics
Analysis of protected area network connectivity serves to identify protected areas that are not well connected to other protected or natural areas, so a requirement for candidate metrics is that values should be highest when surrounding habitat quality is high and safely accessible. Varying movement cost but not mortality with habitat quality does not yield a SAMC disperser distribution metric that meets this criterion because dispersers accumulate where good quality habitat is adjacent to disturbance (e.g. SAMC R high contrast scenario in Fig. 2) so we omitted SAMC R from further consideration.
Least-cost path distances are more computationally demanding than Euclidean distances, making PARC and resistant kernels more costly and difficult to apply over large landscapes. For example, a naive implementation of PARC (e.g. https:// github. com/ LandS ciTech/ LSTDC onnect/ tree/ parc) for a circular park of 500 km 2 analyzed at a resolution of 1 km 2 assuming mean displacement of 40 km requires construction of a distance matrix with 27 million elements. Increasing the mean displacement to 80 km increases the matrix size to 59.8 million elements. Memory requirements can be reduced by iteratively applying the method to different source locations (Compton et al. 2007;Cushman and Landguth 2012) or groups of them, and it would be possible to improve the useability of an iterative approach with parallelization. Other options for reducing the cost of calculations include reducing number of nodes considered (Van Moorter et al. 2022) or simplifying neighbourhood geometry (Drielsma et al. 2007(Drielsma et al. , 2022. We have not pursued these options because we believe there are good theoretical reasons to prefer SAMC metrics over PARC/resistant kernel approaches, including explicit representation of mortality and movement that does not assume knowledge of least cost paths ( Table 1). Analysis of a small subset of protected areas indicated high correlations between PARC and other metrics (Online Resource Fig S7), and we omitted PARC from further consideration.

Analysis of protected area connectivity
For each displacement distance (2, 5, 10, 20, and 40 km) we calculated each raw metric variant (Table 2), except SAMC R and PARC, using habitat layers with and without the human disturbance footprint (Fig. 1). We summarized the spatial information by dividing the human impact layer by the layer without human impacts to remove the effects of natural features and focus on the effect of human footprint. Ratio metric variant scores indicate the proportion of the score that can be attributed to naturally high cost features (i.e., high elevations, steep slopes, large lakes, rivers with high flow, and glaciers), so low values indicate high human impacts. We calculate the mean ratio map pixel value within protected areas to get a score for each metric variant for each protected area, and examine Spearman rank correlations among ratio metric variant scores. We limited our analyses to protected areas larger than 1 km 2 in area. To better understand differences in the behaviour of metrics, we also mapped percentiles of each metric and the mean and standard deviation of these percentile scores across an example region that includes mountains, plains, naturally fragmented coastline, and human impacts of various intensities (British Columbia and Western Alberta). All R (R Core Team 2021) code and data necessary to run this analysis is available at the GitHub repository for this paper: https:// github. com/ LandS ciTech/ LSTD-Conne ctivi ty-Paper/. Spatial inputs and outputs are available for download at https:// osf. io/ bxvd4/.

Results
When mean displacement is low (2 or 5 km) there are high correlations (> 0.85) among all ratio metric variants except SAMC BH2 (Fig. 3). Resistance to moving into low quality areas means that in some cases BH2/5 metrics assign highest values to high quality habitat adjacent to disturbed areas (Figs. S13 and S14). At a national scale, other metrics largely reflect underlying patterns of human disturbance-more highly developed areas are less intact and connected (Figs. 1, S13 and S14). Correlations among ratio metric variants remain reasonably high (> 0.82) at 10 km mean displacement. At larger spatial scales (20 km or 40 km) differences in the behaviour of intactness, SAMC, and simple distance-weighted sum metrics become more apparent (Fig. 3). This pattern of high correlations among metric variants at small spatial scales holds when we examine correlations within ecoregions (Figs. S17 to S31). Linear features in otherwise intact areas have a greater impact on SAMC metrics, highlighting the sensitivity of SAMC metrics to varying road densities in rural areas (Fig. 4). A closer look at spatial patterns in Western Canada shows more clearly that, even when mean displacement is high (20 km), there is high agreement among metric variants in roadless and highly developed areas (Fig. 5). Areas of highest disagreement are naturally fragmented by water or mountains with moderate amounts of human disturbance (Fig. 5).
Comparison of ratio (without naturally high costs) and raw (with naturally high costs) metric values (Figs. 6 and S8-S16) emphasizes the importance of assumptions about landscape cost or quality. In our quality surface, coastal waters and many portions (steep slopes, high elevations, glaciers) of the Rocky and Columbia mountains are considered low quality (Fig. 1), so coast and mountain areas are not well connected even where there is little impact of human activity (Fig. 6 row 2). Calculating ratios to remove the effects of natural features (Fig. 6 row 3) highlights that much of this low connectivity is not attributable to human activity. From a human footprint perspective, the mountain ranges represent relatively intact north-south passages (Fig. 6 row 3). From the perspective of organisms that cannot travel through high mountains, raw metrics highlight lower elevation valleys of interior British Columbia and the east slope of the Rocky Mountains as potentially important north-south movement corridors (Fig. 6 row 2).

Discussion
Globally, resources are increasingly being used to establish and expand well-connected networks of protected areas to meet global targets. Metrics are required to measure progress towards achieving those targets, including for assessing the connectivity of protected areas. With that in mind, we developed fast, parallel implementations of movingwindow, raster-based connectivity metrics, applied them across Canada using 300 m resolution cost surfaces and compared their results. At relatively small spatial scales (2-5 km mean displacement), correlations among metrics are high, suggesting that any of the metrics we investigated will give similar results and simple metrics will suffice, provided that it is reasonable to assume that movement cost and mortality are either constant or correlated with the number of dispersers from each location (Fig. 7). At larger spatial scales (20-40 km mean displacement) in moderately disturbed regions differences among metrics are more conspicuous, implying that assessing connectivity at these scales in these regions requires more detailed information about disperser distribution, behaviour, and mortality risk (Fig. 7). Together, these results suggest different metrics may be required at different scales and landscape contexts, but in some circumstances simple distance-weighted metrics likely suffice. Locke et al. (2019) have proposed a "Three Global Conditions for Biodiversity Conservation and Sustainable Use" (3Cs) framework that divides the terrestrial world into cities and farms, shared lands, and large wild areas, and recognizes the need for different  (Table 2) for 5149 protected areas across Canada. Black outlines identify groups of correlated metric variants identified by hierarchical clustering (hclust function, complete linkage method; R Core Team 2021). When mean displacement is low (2-5 km; cluster B) correlations among metric variants (except SAMC BH) are high. At larger spatial scales (20-40 km), Intactness (cluster A), SAMC (cluster C), and Exponential/ Uniform (cluster D) metric variants are less correlated with one another. At intermediate spatial scales (5-10 km; cluster E) correlations remain moderately high conservation goals and approaches in these different conditions. They identify conditions using the proportion of area within 100 km 2 hexagons in different land use classes. Our distance-weighted sum and SAMC metrics summarize similar human footprint information in somewhat different ways. The metric that most closely resembles the method used by Locke et al. is a simple uniform kernel with mean displacement of 5 km, which represents a circle of ~ 176 km 2 . Moving-window approaches applied at 300 m resolution give a higher resolution view of the landscape than 100 km 2 hexagons, but at a national scale the patterns identified by all metrics reflect the distinctions among city and farms, shared lands and large wild areas proposed by Locke et al. (2019). If the goal is simply to identify where protected areas fall on this spectrum at a spatial scale comparable to that considered by Locke et al., high correlations among all of our metric variants at 5 km mean displacement imply that the choice of metric variant will not have a big impact on outcomes.
In contrast, the choice of spatial scale or mean displacement does alter outcomes. This result reflects the reality that small areas can be functionally intact for organisms that move small distances, and fragmented from the perspective of organisms that move longer distances. Regions dominated by croplands and cities generally lack areas that can be considered intact at large spatial scales, and also lack animals that rely on large intact areas (Hill et al. 2020). In these regions, analysis at finer spatial scales may be helpful to identify relatively intact fragments used by species that persist in landscapes dominated by human activity. Since intactness depends on dispersal capacity there is no single result that would apply equally for all species, but there are circumstances where we can expect similar estimates of relative connectivity across species with a wide range of dispersal capacities. In particular, areas with high metric values at coarse spatial scales represent roadless regions where large wideranging animals that are sensitive to human disturbance persist (Toews et al. 2018). In these regions, a focus on maintaining high intactness at large spatial Fig. 4 Selected ratio metric variants across Canada, assuming 20 km mean displacement. Ratio metrics measure the proportion of a score that can be attributed to natural features by calculating values for a selected metric variant (Table 2) at each location with and without the human footprint (Fig. 1), and taking the ratio of the two values at each location. Low values indicate high impact of the human footprint scales by limiting industrial development to existing nodes and roads, maintaining large migrations, and establishing large protected areas (Locke et al. 2019) will also maintain intactness at smaller spatial scales. We note that this analysis was conducted at spatial grain (300 m) that is finer than many international human impact or footprint analyses (Kennedy et al. 2019;Locke et al. 2019;Beyer et al. 2020;Venter et al. 2020;Hirsh-Pearson et al. 2022), and relevant for conservation planning, but relatively coarse from the perspective of animals that perceive and move very short distances (Albert et al. 2017;Diniz et al. 2020). We did not investigate sensitivity to spatial grain, but it is reasonable to expect that finer resolution analysis would reveal variation in perceived intactness among animals that move 10's vs 100's of metres, analogous to the effects of varying mean displacement from 1 to 40 km on a coarser landscape that are highlighted in this analysis.
In shared lands with intermediate levels of disturbance, particularly those that are also naturally fragmented, we see the highest disagreement among metric variants. Roads and other linear features act as absorbing barriers in the SAMC model, while in the distance-weighted sum model the impact is only a small loss of habitat. The impact of linear features on the SAMC metrics are highest where those features are surrounded by high quality habitat, which happens most often in moderately disturbed "shared lands" that are also often a high priority for biodiversity conservation. The conservation value of "roadless areas" is widely recognized (Ibisch et al. 2016;Psaralexi et al. 2017;Talty et al. 2020;Dietz et al. 2021), and in many cases we expect intact areas identified by  (Table 2). The method for calculating ratio percentiles is shown in Fig. 6. Ecozones (Agriculture and Agri-Food Canada 2014) are outlined in black (labels in b), and 20 km mean displacement is assumed. In this example region in Western Canada, mean ratio percentiles a are high in regions with little human disturbance (e.g. boreal cordillera), and low where the human footprint is high (e.g. prairies) (b). Disagreement (c) is measured as the standard deviation of ratio percentiles among metric variants. In highly developed regions and regions far from human impact (b), there is more agreement among ratio metric variants (c). Areas of highest disagreement (c) are naturally fragmented by water or mountains (d) with a moderate human footprint (b). Zero quality areas are omitted from this analysis (gray) SAMC metrics will be the same roadless areas identified by simpler metrics such as distance from a road (e.g. Ibisch et al. 2016). However, the SAMC model allows a more nuanced view of the impact of roads, which can include variation in movement behaviour and mortality risk among species and road types, and the cumulative effects of multiple roads and other disturbances.
Comparison of results with and without naturally high costs in British Columbia and Alberta highlights the complexity of landscapes that are naturally fragmented by mountains or water, and sensitivity to assumptions about the impacts of these natural features on movement, mortality, and habitat quality (see also Marrec et al. 2020). It is useful to recognize that mountains can serve as refuges and movement corridors for species that can move through high elevations and across steep slopes, while also acting as barriers for species that cannot, and that both low elevation areas and mountain ranges are important for conservation. Indeed, the Yellowstone to Yukon conservation initiative (Hilty et al. 2020) focusses on the Fig. 6 Ratio percentiles are calculated by taking the ratio of raw metric values with (top row) and without the human footprint ( Fig. 1), then converting to percentiles (bottom row). Results are shown for four metric variants across an example region in Western Canada, assuming mean displacement of 20 km. Percentile maps highlight broad agreement in ranking among raw metric variants (middle row), even where the raw values differ (top row). In areas of most disagreement among metric variants (Fig. 5c) on the south coast and southern interior, ratio SAMC metrics give relatively lower scores (bottom row). Raw SAMC metrics highlight the impact of linear features in otherwise undisturbed areas (top row) eastern, more mountainous portion of the montane and boreal cordillera ecozones (Fig. 5), and emphasizes the importance of mountain ranges for connectivity. Dry low-elevation portions of the Okanagan valley (south-western portion of the boreal cordillera ecozone in Fig. 5) have also been identified as a priority for species-at-risk protection in Canada (https:// www. canada. ca/ en/ servi ces/ envir onment/ wildl ifeplants-speci es/ speci es-risk/ pan-canad ian-appro ach. html# toc). It is encouraging that our national-scale species-agnostic analysis reflects recognized priority regions for conservation such as the Yellowstone to Yukon and the South Okanagan.
At a national scale, we believe our results reflect biologically meaningful variation in intactness among regions. Within regions, we are less confident that simple species-agnostic models assuming all undisturbed areas that are not high elevation, ice, water, or steep slopes are high quality, low cost, and low mortality (Figs. 1 and 7) capture enough of the complexity of landscapes and species responses to them. In general, species-agnostic approaches do not reliably provide for the needs of all species (Wood et al. 2021). Migration and dispersal behaviour is often quite different from habitat selection behaviour within home ranges (e.g. Zeller et al. 2012;Vasudev et al. 2015;Roffler et al. 2016;Abrahms et al. 2017;Keeley et al. 2017;Scharf et al. 2018;Tack et al. 2019;Diniz et al. 2020), so habitat quality, movement cost and mortality risk are not always correlated with one another. Accounting for climate variation and climate change can also alter priorities (Littlefield et al. 2017;Schloss et al. 2021). There is good evidence that some species suffer increased mortality near roads (Trombulak and Frissell 2000;Bastille-Rousseau et al. 2011;Mumma et al. 2018;Fabrizio et al. 2019), and that road avoidance leads to loss of available habitat (Trombulak and Frissell 2000;Polfus et al. 2011;Laurian et al. 2012;Dickie et al. 2020). However, the same roads can also facilitate movement or attract animals to roadside food resources (Roever et al. 2010;Bastille-Rousseau et al. 2011;Mumma Fig. 7 Guide for selecting among SAMC and distance-weighted sum metrics (Tables 1 and 2). Green solid arrows and checkmarks indicate "ok to assume" and orange dashed lines with xmarks indicate "not ok to assume". Start in the upper left corner Dickie et al. 2020). A model assuming all types of roads have the same impact on all species is not adequate to capture the complexity of wildlife responses to linear features (Trombulak and Frissell 2000). We also note that the data currently available for roads associated with resource extraction activities in Canada is neither consistent nor complete (Poley et al. 2022), so results derived from this information should be interpreted with caution.
Given the limitations of species-agnostic approaches and available data, we believe our most useful contribution in this work is not maps per se, but tools and insight that we hope will make rasterbased moving-window methods easier to understand, and easier to apply across large landscapes. These methods allow consideration of heterogeneity between and within protected areas, addressing concerns about simple graph-based metrics such as Prot-Conn (Saura et al. 2017). There is a set of simplifying assumptions associated with each available metric, and none are appropriate for all species in all circumstances, but some are easier to interpret biologically than others (Table 1, Fig. 7) (Fahrig et al. 2021;Diniz et al. 2020). Resistant kernel approaches (Compton et al. 2007;Cushman and Landguth 2012;Diniz et al. 2020) such as PARC (Drielsma et al. 2007; Biodiversity Indicators Partnership 2019) include a dissipation function that is somewhat analogous to mortality (Table 1, Eq. 1 and Fig. 2) (Diniz et al. 2020), but difficult to interpret because it is implicit and confounded with movement behaviour. The metrics also assume that dispersers have perfect knowledge of the local landscape that allows them to choose least cost paths while at the same time lacking the ability to discriminate between poor and good quality locations when settling. To understand this point, consider dispersers from a source pixel j. The proportion of those dispersers that end up at a focal pixel i depends on the least-cost distance between i and j, but not specifically on the cost at pixel i. A nearer high-cost pixel will receive more dispersers than a slightly further low-cost pixel. It is difficult to imagine any organism moving and settling in this way. In contrast, the SAMC model (Fletcher et al. 2019) assumes that at each timestep dispersers choose among immediately adjacent pixels, informed by local knowledge of relative cost of movement, and no knowledge of the wider landscape. This assumption is more likely to be reasonable for naive dispersers seeking new habitat or territories than for movement within home ranges (Zeller et al. 2012;Roffler et al. 2016;Abrahms et al. 2017;Keeley et al. 2017;Scharf et al. 2018;Diniz et al. 2020) or along migratory routes (Tack et al. 2019;Sawyer et al. 2020;Fullman et al. 2021). However even naive dispersers may maintain some degree of movement direction over time (Avgar et al. 2016) including passive-wind dispersed organisms (van Putten et al. 2012;Kling and Ackerly 2021). An interesting open question is how results from SAMC models parameterized to reflect species-specific habitat selection, movement preferences and mortality risks compare to agent-based movement models parameterized from integrated step selection functions that include correlations among steps (e.g. Signer et al. 2017;Fryxell et al. 2020;Whittington et al. 2022). We remain intrigued by the possibility of SAMC models serving as a useful compromise between biological realism, simplicity, and computational cost.
There are a variety of innovative and informative methods for estimating connectivity model parameters from telemetry (e.g. Signer et al. 2017;Fryxell et al. 2020;Whittington et al. 2022) and genetic (e.g. Galpern and Manseau 2013;Jasper et al. 2022;Martin et al. 2023) data. However, collecting data to parameterize species-specific movement or connectivity models is generally expensive and timeconsuming, and cannot be done for all species in all contexts. We hope we this work serves to demonstrate how sensitivity analysis can help clarify where additional information might influence decisions (Nicholson and Possingham 2007). In particular, clarifying relationships among distance-weighted sum and SAMC metrics may provide a useful perspective on questions about when and where it is important to consider effects of fragmentation per se in conservation planning (Fahrig 2017;Fletcher et al. 2018;Fahrig et al. 2019). Simple distance-weighted sum metrics (Miguet et al. 2017;Beyer et al. 2020) primarily measure the amount of habitat in the vicinity of each location on a landscape. Metrics constructed from SAMC unconditional short-term disperser distributions (Fletcher et al. 2019) account for habitat amount in a similar way, while allowing for additional effects of the landscape on movement behaviour and mortality. Where distance-weighted sum and SAMC metrics are highly correlated with one another it may not be necessary to precisely characterize effects of the landscape on movement behaviour and mortality when identifying areas of value and concern, and it is likely also difficult to empirically disentangle effects of configuration and habitat amount. Where there is less agreement among distance-weighted sum and SAMC metrics it may be more possible and important to consider impacts of landscape configuration, and additional investment in species-specific data collection and analysis may be warranted.
We did not include Circuitscape and Omniscape methods (McRae et al. 2008(McRae et al. , 2016 in this analysis because our goal was to compare metrics that distinguish intact or well-connected areas from those that are not. Circuit methods yield highest values where movement is constrained and so are more useful for highlighting potential corridors and pinch points than intact areas. We also omitted SAMC R metrics from consideration for this reason. Pither et al. 2021 recently completed an omnidirectional analysis for Canada using the same underlying cost information used in this analysis. An interesting follow-up to this work would be a comparison of SAMC metrics intended to identify areas of high mortality and/or use by dispersers to omnidirectional results and empirical movement, roadkill, and road survey data (Boyle et al. 2017;Fabrizio et al. 2019;Pither et al. 2021;de Rivera et al. 2022).
Patch-free raster-based moving window connectivity analysis methods are computationally intensive which has limited their use and utility thus far. Recent efforts to develop computationally efficient implementations of Circuitscape, Omniscape, and Randomized Shortest Paths in Julia have improved the useability of those methods (Hall et al. 2021;Van Moorter et al. 2022). Conceptually, SAMC metrics provide an intriguing alternative to circuit and resistant kernel metrics, but an implementation that requires construction of a full transition matrix (Marx et al. 2020) is memory intensive and difficult or impossible to apply over large landscapes. Development of an alternative parallel iterative moving-window implementation is a useful step toward improving the useability of SAMC metrics. We recognize there is more work to do on documenting and testing our implementation, and we are exploring options for integrating parallel processing into the samc package (Marx et al. 2020) rather than continuing development of an alternative tool. We hope these efforts will help make SAMC metrics more widely accessible and easily useable, and thereby increase the diversity of landscape perspectives considered in conservation planning.