Evaluating the suitability of close‐kin mark‐recapture as a demographic modelling tool for a critically endangered elasmobranch population

Abstract Estimating the demographic parameters of contemporary populations is essential to the success of elasmobranch conservation programmes, and to understanding their recent evolutionary history. For benthic elasmobranchs such as skates, traditional fisheries‐independent approaches are often unsuitable as the data may be subject to various sources of bias, whilst low recapture rates can render mark‐recapture programmes ineffectual. Close‐kin mark‐recapture (CKMR), a novel demographic modelling approach based on the genetic identification of close relatives within a sample, represents a promising alternative approach as it does not require physical recaptures. We evaluated the suitability of CKMR as a demographic modelling tool for the critically endangered blue skate (Dipturus batis) in the Celtic Sea using samples collected during fisheries‐dependent trammel‐net surveys that ran from 2011 to 2017. We identified three full‐sibling and 16 half‐sibling pairs among 662 skates, which were genotyped across 6291 genome‐wide single nucleotide polymorphisms, 15 of which were cross‐cohort half‐sibling pairs that were included in a CKMR model. Despite limitations owing to a lack of validated life‐history trait parameters for the species, we produced the first estimates of adult breeding abundance, population growth rate, and annual adult survival rate for D. batis in the Celtic Sea. The results were compared to estimates of genetic diversity, effective population size (N e ), and to catch per unit effort estimates from the trammel‐net survey. Although each method was characterized by wide uncertainty bounds, together they suggested a stable population size across the time‐series. Recommendations for the implementation of CKMR as a conservation tool for data‐limited elasmobranchs are discussed. In addition, the spatio‐temporal distribution of the 19 sibling pairs revealed a pattern of site fidelity in D. batis, and supported field observations suggesting an area of critical habitat that could qualify for protection might occur near the Isles of Scilly.

ling approach based on the genetic identification of close relatives within a sample, represents a promising alternative approach as it does not require physical recaptures. We evaluated the suitability of CKMR as a demographic modelling tool for the critically endangered blue skate (Dipturus batis) in the Celtic Sea using samples collected during fisheries-dependent trammel-net surveys that ran from 2011 to 2017.
We identified three full-sibling and 16 half-sibling pairs among 662 skates, which were genotyped across 6291 genome-wide single nucleotide polymorphisms, 15 of which were cross-cohort half-sibling pairs that were included in a CKMR model. Despite limitations owing to a lack of validated life-history trait parameters for the species, we produced the first estimates of adult breeding abundance, population growth rate, and annual adult survival rate for D. batis in the Celtic Sea. The results were compared to estimates of genetic diversity, effective population size (N e ), and to catch per unit effort estimates from the trammel-net survey. Although each method was characterized by wide uncertainty bounds, together they suggested a stable population size across the time-series. Recommendations for the implementation of CKMR as a conservation tool for data-limited elasmobranchs are discussed. In addition, the spatio-temporal distribution of the 19 sibling pairs revealed a pattern of site fidelity in D. batis, and supported field observations suggesting an area of critical habitat that could qualify for protection might occur near the Isles of Scilly.

| INTRODUC TI ON
Estimating contemporary population demographic parameters such as census size (N c ) is fundamental to understanding a species' recent evolutionary history, and is essential for effective conservation of endangered species. The status of marine fish populations is often assessed by estimating their relative abundance based on catch per unit effort (CPUE), obtained from data on commercial landings and/ or fisheries-independent surveys. However, CPUE estimates can be biased by misreporting of catches and variations in catchability of the animals with different sampling gears (Maunder & Piner, 2015).
Data limitations are particularly pertinent to elasmobranchs, assessments of which are often hindered by taxonomic confusion (Iglésias et al., 2010) and insufficient knowledge of their biology (ICES, 2020;IUCN, 2015). These data limitations represent a major conservation concern; many elasmobranchs have suffered population declines and local extinctions as a result of fishing pressure over the past century, and many are still caught as bycatch despite landing bans (ICES, 2020;Simpson & Sims, 2016).
Alternatively, mark-recapture approaches can be used to estimate absolute abundance (Cormack, 1964;Jolly, 1965;Seber, 1965), with the added benefit that they can reveal patterns of animal movement and habitat use. As a result, mark-recapture has been a popular approach to study elasmobranch species (Biais et al., 2017;Corrigan et al., 2018;Neat et al., 2015). However, applications to benthic elasmobranchs such as batoids have generally suffered from low animal recapture and tag recovery rates Bird et al., 2020). For elasmobranch populations facing methodological constraints such as these, novel approaches are urgently needed to estimate recent population trends and inform conservation decisions.
Rapid advances in molecular genetic approaches in the last few decades have provided novel means of assessing population trends. Close-kin mark-recapture (CKMR) has recently emerged as a method of estimating demographic parameters such as absolute abundance, population growth rates, and survival rates (Bravington, Skaug, et al., 2016). CKMR builds on traditional mark-recapture approaches by making use of genetic data obtained from small biopsies (e.g. muscle or fin clips) to identify closely related individuals within a population; here, genotypes can be considered as 'tags' and individuals' relatives in a sample can be considered as 'recaptures' based on the principles of Mendelian inheritance. The principles of mark-recapture are conserved in CKMR, in that a higher proportion of recaptures reflects a smaller population size. Since the method requires only a single capture of individuals, circumventing the need for physical recaptures and the additional stress this would inflict on the animals, it represents a promising approach requiring less time at sea to assess a critically endangered elasmobranch population.
The idea behind CKMR is by no means new. Skaug (2001) and Nielsen et al. (2001) initially proposed that individual genotypes can function as genetic tags from which abundance can be estimated.
However, genotyping and identifying close kin among large populations has only become feasible following recent advances in sequencing technologies, such as next-generation sequencing (NGS), which have enabled the genotyping of large numbers of samples across the many loci (several thousand) required to accurately identify close kin. Thus far, CKMR has been successfully applied to populations of only a handful of species, including southern bluefin tuna Thunnus maccoyii (Bravington, Grewe, et al., 2016), a few salmonids (Prystupa et al., 2021;Ruzzante et al., 2019;Wacker et al., 2021), white shark Carcharodon carcharias , grey nurse shark Carcharias taurus  and thornback ray Raja clavata (Trenkel et al., 2022).
In addition to facilitating a novel demographic modelling approach, the identification of close kin can provide insights on elasmobranch behaviour and habitat use, which can reveal areas that may qualify for additional protective measures. For example, the identification of cross-cohort siblings of speartooth sharks Glyphis glyphis across three Australian river systems revealed patterns of adult movement and breeding behaviour (Feutry et al., 2017), while the high degree of genetic relatedness within aggregations of basking sharks Cetorhinus maximus in the North-East Atlantic revealed kin-associated behaviour and identified potentially important migratory corridors for this species (Lieber et al., 2020).
The blue skate Dipturus batis (Linnaeus, 1758) is a large-bodied rajid with a patchy distribution across the North-East Atlantic Ocean, occurring in higher densities in the Celtic Sea and Rockall (Frost et al., 2020). D. batis was only recently differentiated from its larger congener, the flapper skate D. intermedius, following recent morphological and genetic investigations (Griffiths et al., 2010;Iglésias et al., 2010). Much of the presumed knowledge on both species is being re-assessed and still lacks detail, and as a result, management will continue on the basis of a single common skate complex (D. batis complex) until species-specific assessments improve (ICES, 2020).
At present, both species are classed as Critically Endangered by the IUCN (Dulvy et al., 2006). Current management implies a total landing ban on common skate in EU waters since 2009, though postban landings and bycatch are still reported (ICES, 2020;Simpson & Sims, 2016).
Anecdotal reports indicate high rates of common skate bycatch in the Celtic Sea, an important conservation and socio-economic concern in the context of a multi-species fishery. Initial surveys K E Y W O R D S abundance estimation, blue skate, close-kin mark-recapture, conservation, elasmobranch, fisheries identified that common skate bycatch in the Celtic Sea is dominated by D. batis . In an effort to evaluate current management regimes, longitudinal fisheries-dependent surveys have been initiated in collaboration with British  and French (Barreau et al., 2016) fishing fleets in order to develop demographic models and identify biologically important sites for the species. Initial results from these surveys suggest the occurrence of juveniles and adults in the area, displaying a high degree of site fidelity, while population genomic analysis has identified the occurrence of siblings in the area and limited gene flow with offshore populations (Delaval et al., 2022). Altogether, these findings suggest that the Celtic Sea may host important reproductive and/or nursery sites for D. batis.
The Celtic Sea D. batis population represents a unique opportunity to test the suitability of CKMR as a demographic modelling tool for a data-deficient benthic elasmobranch; the apparently siteattached nature of D. batis in the Celtic Sea, the occurrence of close kin and the longitudinal samples available from standardized surveys make it an ideal candidate for CKMR. In this study, we applied a genome-wide genotyping approach (DArTseq™, Kilian et al., 2012) on samples collected during fishery-dependent surveys in the Celtic Sea between 2011 and 2017 to identify half-siblings and generate the first estimates of the number of breeding adults (N), adult population growth rate and adult survival rate for the population using CKMR. To evaluate the suitability of CKMR as applied to this population, these estimates were compared to molecular estimates of N e (the evolutionary analogue of N, see Waples et al., 2018) and genetic diversity, and to relative abundance estimates (CPUE) obtained from the survey. In addition, the spatio-temporal distribution of close kin would allow us to evaluate skate movements and identify potential areas of biological interest that may qualify for additional protection.

| Sample collection and genotyping
Fishery-dependent common skate surveys were performed in the Celtic Sea in collaboration with the fishing industry in 2011, and from 2014 to 2017. The surveys, which ran in early autumn (August-October), sampled using fixed trammel nets along a transect of stations running 12-80 NM to the south and west of Newlyn, Cornwall, UK ( Figure 1). Additional exploratory stations outside the transect area were also surveyed to assess the extent of D. batis' distribution. Further details on the sampling protocol are described in Bendall et al. (2018). Across survey years, biopsies (fin or muscle clips) were taken from 1140 individuals and stored in 96% ethanol or RNAlater®. Skates were also sexed and measured (total length, cm).
A pilot study was performed to determine the power of relationship inference using DArTseq™ genotyping, and the probability of obtaining parent-offspring or half-sibling pairs among the samples (Appendix S1). The results revealed that parent-offspring pairs were unlikely to be found, as samples mostly comprised young adults, whereas half-sibling pairs (HSPs) could be identified with high levels of precision. Based on these results, we opted for a HSP CKMR approach.
Due to cost considerations, we were unable to genotype all samples. Of those available, 683 were selected for genotyping, of which F I G U R E 1 Sampling locations of all Dipturus batis in this study (black points), full-siblings (red squares), and half-sibling pairs (orange and yellow triangles). Straight lines are drawn between fullsibling (red) and half-sibling (orange) pair capture locations. Half-siblings captured in the same haul are indicated by yellow triangles. Latitude and longitude are in decimal degrees. 387 had already been genotyped by Delaval et al. (2022) using the same method. Samples primarily consisted of juveniles and young adults, but spanned as wide a size range as possible to maximize the number of cohorts included in the model. Genomic DNA was extracted using a DNeasy® Blood & Tissue kit (Qiagen), quantified on a Qubit fluorometer (Thermo Fisher Scientific) and adjusted to 10-60 ng/μl prior to sequencing. To assess the suitability of the DNA for restriction enzyme digestion, we performed a mock sample digest in CutSmart® Buffer (New England Biolabs) for 2 h at 37°C, and resolved all samples on 0.8% TAE electrophoresis gels. In addition, 33 samples were visualized on a Genomic DNA ScreenTape® for a more detailed visualization of DNA quality. DNA was then sent to Diversity Arrays Technology (DArT) Pty. Ltd. for genotyping using DArTseq™ technology.
Genotyping was performed following standard protocols as described in Kilian et al. (2012). DArTseq™ combines complexity reduction methods and NGS platforms, and is optimized for each organism.

| SNP filtering
Single nucleotide polymorphisms were filtered based on a call rate of 95%, and when duplicate loci were present, only that locus with the highest call rate was retained. After this step, the proportion of scored loci per sample was assessed; all samples had a sufficiently high score rate (≥90%) to be retained. Monomorphic loci and those with low minor allele frequencies (MAF < 0.05) were identified using adegenet (v 2.1.2, Jombart, 2008;Jombart & Ahmed, 2011), as implemented in R (v 3.6.2, R Core Team, 2019), and subsequently removed. Next, we tested for conformation of loci to Hardy-Weinberg proportions using the R package pegas Paradis, 2010), performing an exact test based on Monte Carlo permutation of alleles (Guo & Thompson, 1992) with 1000 replicates. After applying the false discovery rate correction method of Benjamini and Hochberg (1995), loci were removed if they deviated significantly from Hardy-Weinberg equilibrium (at a significance threshold of α = 0.05). We then tested for linkage disequilibrium (LD) among loci using the R package snpStats (v 1.36.0, Clayton, 2020) and removed one locus from each pair of loci for which R 2 > 0.80. Because human error could lead to sampling an individual multiple times or to contamination during molecular laboratory work, we looked for duplicate samples based on a threshold of 629 mismatching loci (roughly 10% of remaining loci) using the R package CKMRsim (Anderson, https://doi. org/10.5281/zenodo.3519358). Where duplicates were found (i.e. >90% genetically identical), only the sample with the highest score rate was retained. Following these filtering steps, summarized in Table S1, the resulting dataset contained 662 individuals genotyped at 6291 loci (Table 1, Figure 1). relatedness-finding methods, we also searched for related pairs using ML-relate (Kalinowski et al., 2006) and calculated pairwise relatedness (r) using the Wang estimator in the R package related (v. 1, Pew et al., 2015).

| Close-kin mark-recapture
Close-kin mark-recapture builds upon traditional mark-recapture (MR) abundance estimation by integrating information on the relatedness of individuals. In short, CKMR consists of determining the relationship (in this case as half-siblings or not) of all pairs of individuals in a sample and comparing these to the prior probability of relatedness, given the life-history parameters. Calculating these prior probabilities amounts to building a population dynamics model.
The life-history parameters, if known for the species, can include sex-specific growth rates, age at maturity, and fecundity-at-age. This component of the model is flexible and can be adapted to the species in question by the addition or removal of life-history parameters (Bravington, Skaug, et al., 2016). Models for data-deficient species are by necessity simplified, requiring a number of assumptions to be made. After accounting for the observed related and unrelated pairs inferred from the genetic data, the population dynamics model then feeds into a demographic model from which initial population size and population growth rate can be estimated. In a HSP model, adult survival rates can also be estimated using information on the time interval between siblings' birth years, inferred from their estimated ages.
The HSP approach requires knowledge of the age of the individuals at the time of sampling. When sampling non-lethally in the field, this implies estimating an individual's age based on its size. The datadeficient nature of D. batis means little is known of its life-history traits. Through a mark-recapture experiment, Barreau et al. (2016) were able to use a sclerochronological approach using vertebral growth readings to propose an age-at-length model for D. batis, using the von Bertalanffy (1938) growth equation (Equation 1), which is commonly used to model fish size at age. We estimated the age (in years t) of each individual based on their length (in cm L) using this equation, rounding down to the year (Barreau et al., 2016).
We recognize that the parameters of this equation still require validation from further recaptures (Barreau et al., 2016), and that they lack defined uncertainty bounds that should ideally be included in a CKMR model. However, they offer the best method currently available to estimate age-at-length for D. batis. In order to assess the sensitivity of our CKMR model to erroneous age estimates, we repeated the CKMR analysis described below in a (less computationally demanding) maximum-likelihood framework after re-assigning ages to individuals according to their size (see Appendix S5). Owing to a lack of further available life-history data, we opted for a sex-and age-aggregated population dynamics model (Equation 2). That is, we assumed equal growth rates for both males and females, and constant adult survival rates over time. This approach is similar to that applied to other elasmobranchs to date Hillary et al., 2018), which also lacked precise life-history parameters.
In our HSP population dynamics model (Equation 2), the probability that the relationship K between any two individuals i and j is half-sibship (HSP), given their respective birth years (b) and that j was born after i, is dependent on the number of potential parents (i.e. number of breeding adults N) in the year of j's birth and on the adult survival rate (φ) The rationale behind this formula, which is adapted from Hillary et al. (2018), is firstly that the probability that i's mother (M i , say) survives from b i to b j is b j −b i , and secondly, the probability that j is the offspring of M i is 2 ∕ N b j . The latter follows because there are N b j ∕ 2 candidate mothers at time b j . Finally, the factor 2 comes from also taking the father into account.
The fact that j is M i 's offspring with probability 2 ∕ N b j builds on the assumptions that the population is closed. Indeed, the closed population assumption may be valid, given that migration in and out of the Celtic Sea is apparently limited based on genetic (Delaval et al., 2022) and tagging  data. Regarding the issue of differential probability of capture across space and time, the distribution of D. batis in the Celtic Sea is characterized by biological hotspots and spatially variable CPUE . To minimize any potential bias, we selected samples randomly across survey stations and years. where the proportionality symbol ∝ is used because we omit the normalization constant in Bayes' formula. Further, the term P N t = 0 P(r)P( ) is a non-informative prior on the parameters N t = 0 , r and . Finally, is a product across all pairs (i, j) of individuals in the data. For half-siblings the probability P K ij = HSP | N t = 0 , r, is given by Equation (2), while for unrelated individuals it is 1 − P K ij = HSP | N t = 0 , r, . Note that the parameter N t = 0 enters into Equation (2) through Equation (3). Because the set of pairwise relationships K ij are not statistically independent (e.g. individuals i and j occur multiple times) unless a small portion of the population has been sampled, Equation (4) is not the exact posterior distribution, but corresponds instead to the pseudo-likelihood of Bravington, Skaug, et al. (2016).
The Bayesian machinery requires priors to be assigned to all parameters, and in terms of r and , we used the priors to restrict the range of the parameters to ones that are biologically feasible (i.e. between −1 and 1 for r and between 0 and 1 for ). We (1)

| Genetic diversity and effective population size
Population-and locus-wide summary statistics were obtained using GenAlEx (v 6.5, Peakall & Smouse, 2006. We estimated effective population sizes (N e ) using the LD estimator (Hill, 1981;Waples, 2006;Waples & Do, 2010) in NeEstimator (v 2.1, Do et al., 2014). N e is a theoretical estimator of population size that reflects the degree of genetic drift and, thereby, the evolutionary potential of wild populations; conservation thresholds of N e have generally been set to 500 or 5000 individuals to mitigate the loss of genetic diversity and inbreeding depression, though these thresholds are debated (reviewed in Allendorf et al., 2013). The LD estimate assumed random mating, and we report the results when setting a critical value (i.e. MAF at which alleles should be excluded) of 0.05. Confidence intervals were obtained using the Jackknifeover-individuals method. N e was calculated for the samples overall, and for each sampling year.

| Catch per unit effort
In order to contextualize our estimates of population size and growth trends obtained from CKMR, we estimated relative total abundance (i.e. reflecting abundance of all individuals in the population, not just the breeders) by calculating CPUE, a metric that is more familiar to fisheries management.
We calculated annual CPUE rates from D. batis captured during fishery-dependent common skate surveys performed in the Celtic

| Identification of kin-pairs
Using CKMRsim, we identified three full-sibling pairs and 16 HSPs among the 662 genotyped individuals. When comparing the three kin-finding approaches, CKMRsim, ML-relate, and related, the results were largely consistent. However, ML-relate was less conservative, identifying eight additional HSPs. On assessing the relatedness values (r) among these additional pairs of samples, we noticed the additional related individuals from ML-relate generally had lower r values than the rest (Table S3). By including tests for first-cousins in CKMRsim, we identified 196 pairs of individuals that were likely to be distant relatives (e.g. third-order relatives such as cousins). As

ML-relate does not test for third-order relatives, it may have falsely
identified a number of third-order relatives as half-siblings. A conservative approach was adopted, whereby the three full-sibling and 16 HSPs identified using CKMRsim were retained for downstream analysis.
Of the 19 sibling-pairs identified, in four cases the pair of individuals were captured in the same haul. The other sibling pairs were captured between two and 94 km apart ( Figure 1, Table 2). Two individuals were each involved in multiple sibling pairs. The majority of sibling pairs involved at least one individual captured just west of the Isles of Scilly (Figure 1). These results suggest a pattern of site fidelity that might violate the assumptions of CKMR, in that they could lead to spatially biased sampling (Conn et al., 2020). Indeed, comparing the geographic distance between half-siblings with that between all potential pairs suggested that half-siblings were more likely to be found in closer proximity than expected by chance (Kolmogorov-Smirnov test, D = 0.40, p = 0.01; Figure S5). However, given that our sampling protocol was systematic rather than opportunistic, and covered an area of the Celtic Sea believed to capture most of the range of the population, we considered any spatial sampling bias to be negligible (Conn et al., 2020).

| Close-kin mark-recapture
We identified 16 HSPs that could be used for the CKMR model, which was fewer than the 21 we expected based on the results of a pilot study (Appendix S1). The age of the samples, estimated We underscore that, due to the preliminary age-at-length key for the F I G U R E 2 Mean (red line) and 95% credible intervals (black lines) adult breeding abundance of Dipturus batis in the Celtic Sea, estimated using CKMR in a Bayesian MCMC framework. 100 random iterations from the model are shown (grey lines). Estimates for the modelled cohorts (solid lines) and years following the last cohort in the model (dotted lines) are shown. Note that the mean tends upwards towards the beginning and end of the time-series, which is an artefact occurring from taking averages. Note: Sample sizes (N), observed heterozygosity (H o ), expected heterozygosity (H e ), and fixation index (F) are shown together with their standard error, calculated in GenAlEx (v 6.5, Peakall & Smouse, 2006. Effective population sizes (N e ) and their 95% confidence intervals are also shown, calculated using the linkage-disequilibrium method in NeEstimator (v 2.1, Do et al., 2014) at a critical value of 0.05.

TA B L E 3
Mean genomic summary statistics for Dipturus batis for each sampling year, and overall across years species, age uncertainty could not be meaningfully incorporated into the model, so the uncertainty bounds were likely underestimated.
Testing for the sensitivity of the model to erroneous age estimates, we found that re-assigning ages to individuals had a negligible effect on parameter estimates, but slightly widened the uncertainty bounds ( Figure S4).

| DISCUSS ION
In this study, we evaluated the suitability of CKMR as a demographic modelling tool for a data-deficient and critically endangered benthic Estimates of N e were in the order of ~15,000 individuals, which is above the generally accepted conservation thresholds (e.g. of 500 or 5000, see Allendorf et al., 2013) et al., 2018), and the high-throughput sequencing data from which they were derived may be subject to bias due to violations of the assumption of unlinked loci (Waples et al., 2016).
The CKMR estimates suggested an annual adult survival rate of  (Neat et al., 2015), although their estimate included juveniles which likely experience higher natural mortality rates than adults on account of their smaller size. The higher mortality of D. intermedius might also reflect a greater susceptibility to fishing mortality for this species (Brander, 1981). In comparison with existing elasmobranch CKMR studies, the adult survival rate estimates for D. batis were lower than for adult white shark Carcharodon carcharias  and grey nurse shark Carcharias taurus , both estimated at over 0.90. Such a high survival rate might be expected for large oceanic predators given their higher trophic position. It is worth noting that mortality rates for D. batis discards (accidental catches) in commercial tangle-net and trawl fisheries have been estimated at 38.5% (Ellis et al., 2016) and 33.5% (Barreau et al., 2016), respectively, while skates in general experience mortality rates as high as 45% following capture in UK trawl fisheries (Enever et al., 2009). The potential rates of mortality that can be inflicted on skate populations by commercial fisheries may therefore be consequential, supporting the need for rigorous monitoring and implementation of targeted conservation actions.
There is growing evidence of philopatric, site-attached (e.g. Feutry et al., 2017;Neat et al., 2015), and aggregating behaviour (Lieber et al., 2020;Thorburn et al., 2018) exhibited by elasmobranchs, which may relate to preferential feeding, reproductive or nursery grounds that could qualify as conservation areas. The fulland half-sibling pairs identified in this study were sampled in closer proximity to one another (<100 km apart) than would be expected by chance, even when collected several years apart, indicative of site fidelity. Of particular interest was the high density of sibling pairs that occurred in the waters west of the Isles of Scilly. This location corresponded to an area thought to be 'biologically important' by Bendall et al. (2018) where there were increased catches of immature and mature males and females, and the occurrence of sexually active males and egg-bearing females had been observed. Our results therefore support earlier observations, and identify a site warranting further investigation, one that could perhaps benefit from targeted conservation actions.
This study allows us to evaluate sampling-design consider- Appendix S1).
Furthermore, the site fidelity of D. batis suggests that sampling effort should be increased longitudinally and across several stations in order to avoid sampling litter mates (i.e. siblings born in the same year), which complicate CKMR as they violate the assumption of independent samples (Bravington, Skaug, et al., 2016). In practice, sampling litter mates that are the offspring of an adult with especially high reproductive success in a particular year may downwardly bias the abundance estimate; site fidelity would also increase the probability of sampling a litter mate, especially if these two siblings happen to share a trait conferring increased survival, and therefore capture. This is why we opted for a conservative approach by omitting same-cohort pairs. In contrast, it may be argued that given the level of uncertainty in our age estimates, the same-cohort pairs could have been retained.
Despite the preliminary nature of the age-at-length key used in this study (Barreau et al., 2016), it provides the best available age estimator for D. batis at the time of writing. Accurate age estimation remains a challenge in elasmobranch population assessment, often requiring lethal sampling to obtain vertebral growth readings.
However, novel statistical approaches now provide an opportunity to model age based on individual size, though these require considerable sampling efforts (Régnier et al., 2021).
The variable catchability of elasmobranchs by different sampling gears is also an issue to consider. The trammel-net surveys captured mostly sub-adults and young adults (length range = 66-148 cm; Table 1). The relative lack of small juveniles may have limited the number of cohorts sampled and might explain the absence of parentoffspring pairs in our samples, which would have enabled total adult abundance estimation using a parent-offspring pair CKMR approach (Bravington, Grewe, et al., 2016;Bravington, Skaug, et al., 2016).
In contrast, the bottom trawl data used in the CPUE estimates of Barreau et al. (2016) contained a wide range of size classes (20-140 cm), including an abundance of juveniles. These patterns suggest that alternative sampling methods such as trawling may yield a wider size (and age) range of skates that could be advantageous for CKMR abundance estimation, although they may be more destructive.
In summary, our CKMR and N e estimates provide a first ap- The long-term tagging of D. batis in the Celtic Sea, initiated by Bendall et al. (2018) and Barreau et al. (2016), may provide a useful dataset to validate our CKMR estimates once a sufficient number of individuals are recaptured. Thereafter, more precise demographic models could be developed and applied to the population, and to other elasmobranchs sharing similar life-history traits as D. batis.
Despite the limitations discussed, CKMR has proven to be a relatively cost-effective method of estimating important population demographic parameters for a critically endangered data-deficient elasmobranch; the increasing affordability of sequencing technologies is making it possible to genotype large numbers of samples from single sampling events, while offsetting the increasing costs of the repeated surveys at sea needed to physically recapture individuals (which would also incur additional stress on the animals). In addition to the results of a recent seascape genomics study (Delaval et al., 2022), the present findings improve the status of knowledge for what is likely one of the few remaining large populations of D. batis. Together, these studies demonstrate the additive value of genomic SNP-based approaches to elasmobranch research, which can be implemented in different frameworks (e.g. evolutionary, demographic, and environmental) to address longstanding questions of immediate relevance to conservation.

ACK N OWLED G EM ENTS
We are grateful to the various stakeholders involved in the Celtic

CO N FLI C T O F I NTE R E S T
The authors declare no conflicts of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are openly available in the Dryad Digital Repository https://doi.org/10.5061/dryad. n2z34 tn0g.