Spatial structure of boreal woodland caribou populations in northwest Canada

Local population units (LPUs) were delineated in Canada’s recovery strategy for threatened boreal woodland caribou (Rangifer tarandus caribou). Population viability analyses central to contemporary integrated risk assessments of LPUs implicitly assume geographic closure. Several LPUs in northwest Canada, however, were in part delineated by geopolitical boundaries and/or included large areas in the absence of evidence of more finely resolved population spatial structure. We pooled >1.2 million locations from >1200 GPS or VHF-collared caribou from northeast British Columbia, northwest Alberta and southwestern Northwest Territories. Bayesian cluster analysis generated 10 alternative candidate LPUs based on a spatial cluster graph of the extent of pairwise co-occurrence of collared caribou. Up to four groups may be artifacts in as yet under-sampled areas. Four were mapped LPUs that were conserved (Prophet, Parker, Chinchaga and Red Earth). One small group between Parker and Snake-Sahtaneh known locally as the “Fort Nelson core,” and outside any mapped LPU, was also conserved. Finally, one large group, at >136000 km2, spanned all three jurisdictions and subsumed all of six delineated LPUs (Maxhamish, Snake-Sahtaneh, Calendar, Bistcho, Yates, Caribou Mountains) and part of southern Northwest Territories. These results suggest less geographic closure of LPUs than those currently delineated, but further analyses will be required to better reconcile various sources of knowledge about local population structure in this region.


Introduction
Understanding spatial distributions of organisms and the consequences for conservation policy and management decisions remain important challenges (Gaston, 2003;Gaston & Fuller, 2009;Kowalchuk & Kuhn, 2012). Among species assessed by the Committee on the Status of Endangered Wildlife in Canada (COSEWIC), a majority are decided in consideration of species' spatial distributions, i.e., extent of occurrence and area of occupancy, because data for other attributes are typically less available (e.g., Trout, 2013). For purposes of recovery planning, local population units (LPUs) of the boreal population of woodland caribou, Designatable Unit (DU) 6 of Rangifer tarandus caribou (COSEWIC, 2014), were delineated (Environment Canada, 2011, 2012. Where possible, these were informed by available telemetry data (Environment Canada, 2011).
Understanding species' spatial distributions and population dynamics of associated demographic units, at scales less than a species' entire geographic range, is potentially complicated by movement of organisms among spatial subunits over space and time. Population dynamics of a spatially-defined group of organisms depend on vital rates, i.e., births and deaths occurring within the defined area, and immigration and emigration to and from it. Estimation of population dynamics is simplified if immigration and emigration are negligible, i.e., the population unit is geographically closed. Uncertainties in estimating population dynamics of geographically open populations, for which immigration and emigration rates are not known, may compromise conservation policy and management decisions, but methods to directly incorporate immigration and emigration rates are often not immediately practicable. Alternatively, it may be advantageous to develop a method to cluster organisms into groups for which population dynamics might reasonably be assumed to be largely a function of births and deaths, and less so immigration and emigration. Environment Canada's Science Assessment (Environment Canada, 2011) developed an Integrated Risk Assessment (IRA) protocol to assess self-sustainability for boreal woodland caribou LPUs. The assessment integrated three components for which geographic closure was implicit: (1) the probability of short-term, future population growth; (2) the probability of long-term extirpation; and (3) empirically estimated historical population growth. The IRA was applied to 51 LPUs across the Canadian boreal forest, the boundaries of which had been mapped in a variety of ways. Some delineated LPUs could be assumed to be almost certainly geographically closed. However, other LPUs were delineated as broad areas of caribou occurrence without clarity about the extent to which they might comprise smaller LPUs. Still others defaulted to jurisdictional or administrative boundaries. A federal Recovery Strategy (Environment Canada, 2012) and Action Plan (Environment and Climate Change Canada, 2017) called for development of approaches and standards to improve identification of LPU range boundaries.
Approaches and standards to identify LPU range boundaries should be transparent, repeatable, and, to reliably estimate self-sustainability, generate plausible LPUs that are geographically closed. We first explored familiar clustering methods (e.g., Taylor et al., 2001;Shuter & Rodgers, 2013) using telemetry data from Alberta (AB), British Columbia (BC) and Northwest Territories (NT), but these did not generate plausible, geographically closed groups of caribou . Here, we describe a method for grouping caribou into plausible candidate LPUs that may better approximate geographic closure than the existing LPUs.

Study area
The focus of this project was northwest AB, northeast BC and southwest NT ( Figure 1) where empirical evidence demonstrated that caribou move among the three jurisdictions (Kelly & Cox, 2011;Larter & Allaire, 2015;Government of Alberta, 2017;Culling & Culling, 2017), but where some LPUs were delineated on an interim basis to follow provincial and territorial boundaries (Environment Canada, 2011).

Data assembly
We assembled caribou telemetry data collected between 1982 and 2016 in northwest AB, northeast BC and southwest NT. We screened data for spatial errors and used a Lambert Conformal Conic single parallel projection (centre latitude = 60° and centre longitude = -120°), which preserves relative Euclidean distances between points (Taylor et al., 2001;Shuter & Rodgers, 2013). We mapped the first recorded telemetry location of all caribou and visually examined the distribution to identify potential spatial biases in collaring effort which might affect clustering results (Shuter & Rodgers, 2013).

Clustering method
We developed an innovative method to cluster caribou based on spatial relationships that we refer to as a spatial cluster graph. We used individual caribou as random variables (i.e., col-umns in our data table), spatial grid cells as observations (i.e., rows in our data table), and the presence or absence of a caribou in each grid cell as our variates. This is in contrast to other studies with similar objectives that considered individual animals to be observations and grid coordinates to be random variables (e.g., Taylor et al., 2001;Shuter & Rodgers, 2013). Grid cells were set to 10 km 2 to balance the trade-off between spatial precision (where decreasing cell size is better) with ability to detect co-occurrence of caribou (where increasing cell size is better). Testing for sensitivity of results to variation in size of grid cells revealed little effect at 5 km 2 and little spatial structure at 100 km 2 .
Spatial relationships among caribou were characterized by the extent of mutual information shared by all pairs of caribou included in the sample. Mutual information can be thought of as a non-parametric equivalent of a correlation coefficient that quantifies the amount of information (measured in units of "bits," where 1 bit is needed to represent a variable that can take two values, such as a coin toss) obtained about one variable by knowing the value of another variable (Cover & Thomas, 2006). Formally, mutual information (I) between variables X and Y is defined as where H is the "entropy" for the discrete binary distribution of caribou X among grid cells (and equivalently for caribou Y) and is defined as where P(x) is the proportion of grid cells occupied by caribou X based on its marginal distribution. The unit of measure of H is bits, which is why we use base 2 of the logarithm. Mutual information is the difference between the marginal entropy of caribou X and the entropy of caribou X conditional on the entropy of caribou Y. In other words, the mutual information of caribou X and caribou Y provides a measure of how much our uncertainty in the spatial distribution of caribou X is reduced by knowing the spatial distribution of caribou Y.
We then built a graph representing the spatial relationships among caribou, where each "node" of the graph represented an individual caribou and the network topology (i.e., edges between caribou nodes) represented close spatial relationships among caribou pairs in a way that minimized the Minimum Description Length (MDL) score of the graph (Lam & Bacchus, 1994): MDL is the sum of the bits (description length or DL) required to represent the network graph B and associated probabilities, and the bit required to represent the dataset D given the network B.
Network structure was learned as a maximum weight spanning tree based on the MDL scores of candidate networks iteratively fitted from different pairwise combinations of caribou. The resulting network maximized the mutual information among caribou pairs (the weights in the spanning tree) while excluding weakly connected caribou from the network based on the MDL score.
Caribou were then merged together iteratively into groups based on a derivative of mutual information, formally, Kullback-Leibler Divergence, which compares the difference between the probability distributions for pairs of caribou (Polani, 2013). The minimum number of groups occurred where caribou no longer effectively shared mutual information (i.e., complete spatial separation, I = 0), but any number of groups could be defined by setting a minimum divergence between caribou pairs to be considered in the same group. As a result, we were also able to explore solutions with larger numbers of groups. All analyses were completed using Bayesi-aLab 6 (Bayesia S.A.S., Laval, France).

Results
The telemetry data set included over 1.2 million locations collected from 1226 radio-collared caribou (569 GPS and 657 VHF) between January 1982 and November 2016 (Table 1). Adult females are typically the age-sex class radio-collared, and we found no evidence in the datasets that males were tracked (i.e., only females were identified where sex was noted). Minimum convex polygon home ranges averaged 875 km 2 (range 13-5,563 km 2 ) among the 355 caribou with>1,000 GPS or >100 VHF locations recorded. Based on the distribution of first locations, we considered sampling to have been relatively uniform throughout the project area (Figure 2), despite acknowledged gaps that resulted from spatial and temporal variation in collaring effort due to practical constraints and/or to respect concerns among First Nations about collaring.
The spatial cluster graph included 1124 edges among caribou pairs (Figure 3). Ninety-two caribou were excluded on the basis of low MDL scores; they were associated with few telemetry locations and/or small and isolated home ranges. The minimum number of stable groups of caribou identified by the spatial cluster graph approach was 10 and these groups overlapped little spatially (Figure 4). Up to four groups may be artefacts in as yet under-sampled areas. Four were delineated LPUs that were conserved (Prophet, Parker, Chinchaga and Red Earth). A small group between Parker and Snake-Sahtaneh known locally as the "Fort Nelson core," that occurs outside of any mapped LPU, was also conserved in the analysis. Finally, one large group, at >136000 km 2 , spanned all three jurisdictions and subsumed all of six delineated LPUs (Maxhamish, Snake-Sahtaneh, Calendar, Bistcho, Yates, Caribou Mountains) and part of southern Northwest Territories). Arbitrarily increasing the number of groups to 15, for example, revealed finer-scale structuring within the 10 groups ( Figure 5), but reduced confidence that the groups were geographically closed. Specifically, the Chinchaga LPU split into two groups, and overlapping groups were resolved within the Calendar, Bistcho, Caribou Mountains and southern Northwest Territories LPUs.

Discussion
Entropy based clustering methods have a long history in biology (e.g., Dehmer & Mowshowitz, 2011); however, this is the first application of this approach to the problem of describing distributions of wide-ranging ungulates. In contrast to previous cluster analyses that used median locations (Taylor et al., 2001;Shuter & Rodgers, 2013;Wilson et al. 2017), the spatial cluster graph approach permitted use of the en-  3. Spatial cluster graph for boreal woodland caribou (n = 1226), illustrating the structure of the maximumweight spanning tree and the resulting 10 groups. Detail of the connections for the west Great Slave Lake group (northeast of Fort Providence) is illustrated in the inset. Each node (circle) is an individual caribou. Edges (links) between nodes denote pairs of caribou sharing the most mutual information, with thicker edges representing closer spatial relationships (i.e., the weights in the spanning tree).
tire sample of telemetry points. Isolated groups of animals on the periphery of the study area were resolved, such as the Parker and Prophet LPUs in British Columbia, correlating well with known caribou movements and distributions in that area. The 10-group spatial cluster graph grouped caribou from several LPUs into one large range, extending from north of the Mackenzie River to the Snake-Sahtaneh LPU in the south, and east to the Caribou Mountains. This range did not break along portions of some major river and/or road corridors, such as the Mackenzie and Hay Rivers, but did along the Liard, Fort Nelson and Fontas Rivers (Figure 4). Resolving more social groups revealed finer spatial structuring, but some groups overlapped extensively and neither did they break along some river corridors as did others.
We offer that these results contribute improved information about the extent of geographic closure important to the assessment of self-sustainability according to Environment Canada's IRA (Environment Canada, 2011). Nevertheless, some caveats are warranted with regard to both the method itself as well as the Fig. 4. Results of cluster analyses based on the spatial cluster graphs for the 10-group solution. Each colour represents a distinct group. Points noted in the text are labelled: 1) Mackenzie River; 2) Hay River; 3) Liard River; 4) Fort Nelson River; 5) Fontas River; 6) Great Slave Lake; 7) Trout Lake; and 8) Edéhzhíe. results. We group these caveats with respect to how they might bias the results to Type I or II errors, i.e., resolving groups of caribou that are not actually distinct, or not resolving groups of caribou that actually should be distinct, respectively.
With respect to Type I error, limited sampling in some cases may have contributed to resolving distinct caribou groups where there actually may not be, such as the areas immediately south and west of Great Slave Lake (Figure 4), where sampling only began in 2015 and is unlikely to adequately represent the full range of caribou distribution in those areas. In areas east of Trout Lake and on the Edéhzhíe in the Northwest Territories (Figure 4), practical constraints and/or respect for First Nations' concerns limited collaring effort. To the south there is perhaps greater confidence that caribou LPUs may be relatively more isolated and closed, as reflected by the analysis, but even there, Cree elders reported use of habitats by caribou between, for example, the southeast portion of the Caribou Mountains LPU and the Red Earth LPU (Schramm & Krogman, 2001;J. Webb, pers. comm.). The spatial cluster graph resolved the Red Earth and Caribou Mountains LPUs in the absence of telemetry data from between these areas. Similarly, linking consecutive telemetry location points has been unable to demonstrate caribou movements between Red Earth and Caribou Mountains LPUs (Government of Alberta, 2017). Further work to update LPUs is warranted following efforts to obtain more telemetry data from under-sampled areas that are potentially occupied by caribou. Regardless, it will remain important to understand also to what extent different sources of information about species' LPUs at alternative geographic scales, such as that from population genetics and Indigenous knowledge (e.g., Priadka et al., 2018;Schramm & Krogman, 2001) might be reconciled with telemetry-based information.
Several factors potentially affect Type II errors (i.e., not resolving groups of caribou where they should have been), including spatial errors in telemetry data, infrequent long-distance movements, age and temporal resolution of data, and narrow barriers to movement that bisect spatial grid cells. We address these factors in turn. First, the joint probability distribution on which the analysis was based is binary, so groups may be sensitive to spatial errors and/ or occasional long-distance movements (i.e., forays) by individual caribou. Occasional longdistance movements are typically insignificant with regard to identifying LPUs for purposes of management decisions to recover caribou in the nearer term. They play a more important role in the longer term in maintaining gene flow among groups of caribou over long time scales. Further, the sampling has been restricted to females, and the potential importance of male intergroup movements to persistence of caribou LPUs is unknown.
Second, some of the telemetry data was old, such that we might have detected fewer caribou groups than was actually the case by including locations from areas caribou may no longer occupy. Our spatial cluster graph analysis was developed to substitute where traditional cluster analysis did not perform satisfactorily. In this context, it is important to use a large set of locations over space and time. The method is distinguished from social network analysis, which also clusters animals based on spatial relationships, but for which it can be more readily assumed that animals interact socially. For example, Peignier et al., (2019) used a social network analysis (Robitaille et al., 2018) to test among factors affecting space use by caribou at fine spatio-temporal scales, i.e., a true social network, where caribou were defined as interacting if they were less than 50 m apart within 5-min intervals. Instead, we contended with multi-generational data characterized by inconsistent monitoring periods and relocation frequencies. Caribou not even alive at the same time could be members of the same group occupying a distinct area. Even so, this is a property also of traditional cluster analyses that are regardless accepted, for example, to resolve groups of polar bears (Ursus maritimus, Taylor et al., 2001). Nevertheless, the dataset was dominated by recent telemetry (Table 1), and groups were resolved along the southern edge of the study area despite that older data from these areas was included in the analysis. Further analysis of temporally resolved telemetry data, which might provide information about the extent to which landscape change affects change in spatial structure of groups of caribou, was beyond the scope of this project and, in any case, restricted by the comparatively low point densities of the older VHF data. Further research into questions of this nature would be warranted as still longer time series of GPS data become available.
Third, related to the point above, temporal resolution of the telemetry data may have resolved groups of caribou in another respect. Wilson et al. (2017), using familiar cluster analysis, temporally partitioned the data, but found little effect on spatial structure. We did not consider further the potential extent of geographic segregation of caribou throughout the year that otherwise co-occur, for example, during breeding. Neither could our analysis contend with the possibility of behavioural segregation, such as assortative mating, among caribou within distinct social groups despite spatial overlap of those groups during breeding, although both situations may be relevant to defining LPUs and should be further investigated.
Fourth, our spatial cluster graph analysis may be insensitive to narrow barriers to caribou movement (e.g., rivers, highways). Telemetry points from different caribou occurring in a single grid cell but otherwise separated by a physical barrier bisecting the cell is considered evidence of spatial coincidence and increases the likelihood that caribou would be grouped together. This may have contributed to a lack of resolution of groups separated by semi-permeable barriers to movement, e.g., the Hay River and portions of the Mackenzie River in the Northwest Territories, where linking consecutive telemetry location points of caribou reveals non-overlapping use by caribou on either side of the river. Traditional knowledge studies have also reported discrete groups of boreal caribou on either side of the Hay River (Dehcho First Nations 2011, Gunn 2009. Results need to be interpreted in relation to these types of barriers. Population genetic studies similarly revealed greater spatial structure among caribou than that currently recognized by Environment Canada (2011, 2012; within the contiguous Northwest Territories LPU, three broad genetic clusters were detected (Manseau et al., 2017).
The potential for Type II error notwithstanding, some traditional use studies report weaker spatial structuring among some groups of caribou than outlined by Environment Canada (2011, 2012. For example, Dene Tha' First Nation elders do not consider the Yates, Caribou Mountains, Bistcho, Snake-Sahtaneh, Calendar, and southern parts of Northwest Territories LPUs to be distinct. Dene Tha' elders indicate that animals move among Tathlina and Trout Lake (NT), the Cameron Hills and east of the Hay River (AB/NT), along the Yates River, and to north of Bistcho Lake (AB), the Shiekielie River (AB/BC) and to Calendar Lake and the Snake-Sahtaneh River (BC); they are not known to mix with the Chinchaga LPU (M. Munson, pers. comm.).
LPUs constitute logical foci for recovery planning for boreal woodland caribou (Environment Canada, 2011). To the extent that LPUs may not be geographically closed to immigration and emigration of caribou among them, population viability analyses based on birth and death rates alone may contribute to inaccurate assessments of risks posed to LPUs.
Our results contribute in particular to a better understanding of the spatial structure of LPUs otherwise bounded geopolitically, and perhaps other LPUs in the region. As such, this type of analysis can help to position landscape-scale planning for boreal woodland caribou on firmer ecological footing regarding distributions of caribou within and among jurisdictions and can also foster and inform inter-jurisdictional collaboration to better assure self-sustaining caribou populations.