Method selection affects the estimates of residency and site fidelity in bottlenose dolphins: testing sensitivity and performance of different methods using mark-resight data

Metrics

R and SF assessments were first separated based on their actual use as expressed by the authors in peer-reviewed literature. However, given the wide variety of criteria among studies (see Huesca-Domínguez et al., 2024), we started with the simplest concept and thus considered SF indices as those using sighting frequencies regardless of temporal scale and/or sampling periodicity to classify individuals (Table 1). The latter resulted in 11 indices to assess SF, based on different types of variables (i.e., categorical in nominal or ordinal scale, or numerical).

Conversely, definitions in which individuals were directly classified using different measures of occurrence within a given timeframe (monthly, temporary, seasonal, and/or annual) were defined as R indices. Thus, eight R indices corresponded to criteria provided by peer-reviewed publications that defined different classifications of individuals (Table 2). All of these contained a categorical scale with two (residents and non-residents) or three (residents, non-residents, and occasional/transient individuals) types of individuals based on their mark-resight histories.

Description of the model

Mark-resight models to estimate population size, survival and other demographic variables (Nichols, 1992) are based on the idea that each member of the population has a sighting history Y_i = (Y_i1, …, Y_it), where Y_it = 1 if individual i, is sighted at sampling occasion t, t = 1, …, T, and otherwise Y_it = 0 (see Morteo, Rocha-Olivares & Morteo, 2012). One of the most common models is the Jolly-Seber (J-S) (Schwarz, 2001), which has been widely used in a broad range of taxa (e.g., Vargas-Fonseca et al., 2020; Bolaños Jiménez et al., 2021; Ali & Rajamani, 2023; Blázquez et al., 2023), and thus developed, parametrized, extended, and implemented by different researchers across years and study areas.

We chose the algorithm developed by Kéry & Schaub (2012) in the R-project language (simul.js) since it generates sighting histories under the implemented assumptions of the J-S model. The model is based on parameterizations of superpopulations (i.e., the number of animals that use the study area at some point during the study) with constant sighting and survival probabilities and time-dependent entry probabilities (see Williams, Frederick & Nichols, 2011; Bolaños Jiménez et al., 2021).

Following the mark-resight methods, the algorithm simulates sampling strategies based on the identification of individuals (marking) in a target population, to take subsequent samples under the same conditions and determines what fraction of these were identified on different occasions (resightings) (Seber & Schofield, 2019).

The probability that an individual that has not been present since the previous occasion enters the population on occasion t is b_t(t = 1, …, T) and is called the entry probability (Schwarz & Arnason, 1996). Entry could result from on-site recruitment (individuals born locally) or immigration (individuals born outside the study area entering at some point). To ensure that all individuals (N) enter the population at some point during the study, the sum of all entry probabilities must be equal to 1 (Kéry & Schaub, 2012).

The algorithm provides simulated sighting histories of the N individuals, along with the number of individuals that entered the population for the first time and the actual size of the population on each occasion, according to the input parameters provided by the user.

Calibration by simulating scenarios

Calibration was obtained by identifying the optimal values of the parameters from the model, such that its fit was maximized with respect to a set of simulated data (Campbell, 2006; Trucano et al., 2006). In this sense, we used fixed and variable values for parameters to build sighting histories of individuals in simulated wild dolphin populations (see Morteo, Rocha-Olivares & Morteo, 2012). Input parameters for the algorithm (simul.js) were: (1) size of the superpopulation (N), (2) probability of entry (b_t), (3) resighting probability (p), (4) apparent survival probability (Phi), (5) number of marked individuals (n), (6) number of sampling occasions (t_i) and (7) sampling frequency (m). Since dolphin populations around the world can range from very small to over 400 individuals, we used a superpopulation size of 600 individuals (Morteo, Rocha-Olivares & Morteo, 2012; Caruso et al., 2024), which resulted in random sizes of marked populations ranging from a few individuals to over 400, according to the recapture and survival probabilities tested here. The probability of entry (or recruitment) consisted of a random value for each sighting history varying between 0.05 and 0.15; this range was used to represent stochastic variations in the population and sampling effort (Kéry & Schaub, 2012; Morteo, Rocha-Olivares & Morteo, 2012), given that the proportion of the population that is first identified is rather small, but increases over time with an asymptotic trend once sufficient time and sampling effort have occurred (Morteo, Rocha-Olivares & Morteo, 2012).

Resighting probability (p) between sampling occasions was set from 0.01 to 1.00 in increments of 0.10, to represent naturally variable scenarios, considering individuals present only once (i.e., 0.01), those present on different occasions in the study period, and also animals sighted in all samples (i.e., 1.00). On the other hand, survival probabilities (Phi) varied between 0.70 and 1.00 (with 0.05 increments) and were considered to represent the biological effect of sex and age in bottlenose dolphins due to survival heterogeneity between males and females, as well as among developmental stages (i.e., newborn, calf, juvenile, and adult). Survival values reported in the literature are typically higher in females than in males, and lower for calves and juveniles, with trends showing values above 0.7 in all cases (Currey et al., 2009; Ludwig et al., 2021; Arso Civil et al., 2019; Bolaños Jiménez et al., 2021; Barratclough et al., 2024; Schwacke et al., 2024). Also, the sampling duration was set at three years, considering the mean duration obtained in the reviews by Morteo, Rocha-Olivares & Morteo (2012) and ensuring full-time coverage to assess residency using indicators that require at least one or two years of sampling (Huesca-Domínguez et al., 2024). Finally, sightings were clustered using monthly (i.e., 36 occasions), seasonal (i.e., 12 conventional seasons), and annual (i.e., three years) periods to make all calculations.

Sensitivity analysis

In the sensitivity analysis, the input parameters of the model were changed to determine the effect on the behavior of the output values (Sargent, 2010). To this end, sighting histories were originally set to 10,000 simulations following Morteo, Rocha-Olivares & Morteo (2012). The sensitivity of the different definitions of R and SF within the mark-resight model were assessed by changing survival (Phi) and resighting (p) values within the given ranges (0.01–1.00 and 0.70–1.00, respectively) for each unique combination of the fixed values of the model (i.e., population size, entry probability, and number of sampling occasions) to obtain the proposed sets of simulated sighting histories with the same features. Thus, in each of these simulated scenarios, we applied the most common definitions of R and SF to the same data (Tables 1 and 2). We emphasized the variables originally proposed by Ballance (1990), for being among the first and most widely used in the literature (Morteo, Rocha-Olivares & Morteo, 2012; Huesca-Domínguez et al., 2024). Subsequently, we used the estimates of each index for the same simulated data, and averaged the proportion of individuals within each R and SF category to analyze the sensitivity of the estimators to variations of the implemented parameters (i.e., resighting and survival); that is, if the curves of a given index under different values of p and Phi do not change, then the index is insensitive to changes in those parameters over the given range of values. Conversely, if the curve changes drastically as the level of the parameters changes, then the estimator is sensitive to that parameter.

Case studies

We used individual sighting histories from two long-term photographic-identification (Photo-ID) studies on wild populations of Tursiops spp. (Savannah, USA, and Alvarado, Mexico), to compute empirical metrics of residency and site fidelity. Based on the methodology proposed by Morteo, Rocha-Olivares & Morteo (2012), we used data from three consecutive years of study for each population. Sampling efforts from both populations had different frequencies which allowed us to assess the performance of the indices under different empirical scenarios.

The first case study corresponded to data from the common bottlenose dolphin (Tursiops truncatus) photo-ID by the Marine Mammal Laboratory at Universidad Veracruzana (LabMMar-IIB-UV) from May 2006 to August 2010 (Morteo, 2011; Morteo, Rocha-Olivares & Morteo, 2012; Morteo, Rocha-Olivares & Abarca-Arenas, 2014; Morteo, Rocha-Olivares & Abarca-Arenas, 2017). For this research, monthly sighting histories from January 2007 to November 2009 were used, and the study area was covered on each survey, but the dataset contained months during which no sampling was conducted. Nevertheless, information was available for all seasons of the year. In this case, 206 individuals were identified in 29 months across the three years.

The second case corresponded to data from the Tamanend’s bottlenose dolphin (Tursiops erebennus) photo-ID in estuaries from the southern Savannah River to the northern Ossabaw Sound, in Georgia, USA, from April 2009 to July 2017 (Perrtree, Kovacs & Cox, 2014; Kovacs, Perrtree & Cox, 2017). The sampling period covered January 2014 to December 2016, with monthly consecutive surveys. Efforts were made to cover the study area across the surveys, and in total, 455 individuals were identified within the 36-month sampling period, with variable sampling effort, which was greater during the summer months (May–August). However, information was also available for all seasons of the year.

Validation test

Validation aids in determining the degree to which a model is an accurate representation of the empirical wild populations from the model’s intended perspective (Trucano et al., 2006). In this sense, the results obtained in the simulations under the different scenarios were compared with the empirical data of bottlenose dolphin populations from Alvarado, Mexico and Savannah, USA.

First, the individual sighting histories from the empirical data were categorized using the definitions of R and SF previously described (Tables 1 and 2). Simulations with similar resighting and survival values to the empirical data were then used for comparisons. For this, we used the apparent survival (Phi) (i.e., the probability that any of the marked individuals survives and returns to the study area) values for each population as computed under the J-S approach in SOCPROG v2.9 software (Whitehead, 2009). It is noteworthy that both Savannah, USA, and Alvarado, Mexico dolphins showed the same survival (Phi = 0.92) which we used to compare with the simulated estimates of R and SF.

The available data on annual sighting probability for the Alvarado population by Bolaños Jiménez et al. (2021) was also used as a reference (p = 0.88); in this case, the high value refers to individuals having enough resightings such that they may be classified by any of the definitions of R and SF applied in this study, thus it logically rules out transient individuals. The differences between the simulations and the empirical data were then assessed using a contingency table, using the Chi-square or the Mann–Whitney test where appropriate.

Performance

Metrics to assess R and SF utilize different criteria for individual occurrence under different measurement scales (i.e., monthly, seasonal, and annual). Therefore, we used the Gower (1971) distance to compare similarities among individual sighting histories. This measure determines how different two observations are in the presence of mixed data, that is, a combination of numerical and categorical information (Pavoine et al., 2009; D’Urso & Massari, 2019).

The definitions of R and SF considered the different variable types (e.g., numerical and categorical) and measurement scales to calculate the Gower distance. As shown in Tables 1 and 2, the multi-categorical scale was evaluated as a “fuzzy” variable (i.e., several columns were used to define the different possible levels of R and SF in each annual period, such that an individual could belong to different categories throughout the study).

Therefore, once the distance matrix was obtained, the contribution of the result of each index to the general Gower distance was calculated, based on the correlations between the distances obtained for each variable and the global distances obtained by mixing all the variables (Bello et al., 2021). In this sense, the contribution of each index to the dissimilarity was calculated using the “k.dist.cor” function of the “ade4” package (Thioulouse et al., 2018) in the R-project programming language (R Core Team, 2024). These contributions helped to identify indices that have the best resolution or structure for the assessment of residency and site fidelity to determine the similarity among individual sighting histories.

These indices were first applied to data from the two empirical populations (i.e., Alvarado, Mexico and Savannah, USA) and then plotted with the simulated data previously described, for comparison to the simulated scenarios (i.e., when resighting and apparent survival probabilities are both low and/or high) using a violin plot. These results provided the expected statistical distribution of the Gower distances for the selected indices of residency and site fidelity under all the simulated scenarios as a measure of their performance related to the observed (i.e., empirical) values. All statistical analyses were conducted using R-project (R Core Team, 2024).

Calibration

The initial methodology of 10,000 simulations proved computationally intensive given our multiple combinations of parameters for the scenarios, thus the number of simulations was progressively reduced until empirical results showed that output values for the analyzed variables were stable; this occurred at approximately 100 runs, which was set as the optimal number for all the experiments. We also made trials using a range of superpopulation sizes (i.e., N = 100-600) and noted that results did not alter the trends in the data, as these were expressed in proportions; however, confidence intervals narrowed as the number of simulated individuals available for sampling increased; thus, we used the largest N. As expected, the values for the parameters: superpopulation size (N = 600), entry probability (b_t = 0.05 to 0.15), resighting probability (p = 0.01 to 1.00), apparent survival probability (Phi = 0.70 to 1.00), sampling occasion (t_i, I = 1 to 36), sampling frequency (m = once a month) and sampling duration (three years), all conditioned the number of marked individuals (n) available for sighting, and thus affected the computations of R and SF.

Sensitivity

The results for Ballance’s variables are represented in Fig. 1. The number of sighting records (O sensu; Ballance, 1990) showed a linear pattern under the different values of capturability (p), where the higher the survival (Phi), the greater the effect (slope) (Fig. 1A), and none reached the maximum value of 36 sightings (i.e., present in all sampled months). In the case of the permanence, expressed as the number of months between the first and last sightings (P sensu; Ballance, 1990), all curves showed a different asymptote for the selected range of survival probabilities (Fig. 1B). The mean number of months between consecutive sightings (I = inverse of the periodicity, sensu; Morteo, Rocha-Olivares & Morteo, 2012) had slight deviations in the slopes of all scenarios in the simulated data, which yielded in significant differences as the capture probability (p) increased (Fig. 1C).

Only eight of the evaluated indices classified individuals in some degree of R, and these were run in the simulations showing contrasting results. Classifications by Rosel et al. (2011), Dinis et al. (2016), and Conway (2017) only had two types (degrees) of R. For instance, estimates under Rosel et al. (2011) classification (Fig. 2A) showed a similar pattern to the definition by Conway (2017) (Fig. 2C), managing to identify about 25% of resident individuals in unlikely favorable conditions (e.g., p = 1.00, and Phi = 0.95). However, when survival and catchability are at their maximum, about 50% of the individuals may be classified as residents and non-residents. Also, under more probable conditions (i.e., with intermediate values for resighting probabilities) it was possible to distinguish a reasonable proportion of resident individuals (i.e., 25%), but only under high survival probabilities (>0.95).

The simulations under the classification by Dinis et al. (2016) showed that both in poor (i.e., low probabilities of resighting and apparent survival), and in optimal conditions (i.e., where p and Phi are close to 1.00), results were not very different, except for populations with maximum survival and capturability above 0.40 (Fig. 2B). Thus, this index is very strict in classifying resident individuals, managing to include at best 25% of them as residents under the most optimal scenarios.

Simulations for Zolman (2002), showed that the proportion of individuals in the different degrees of R clearly vary under different values of survival. However, individuals were more likely to be classified as transient compared to the other categories despite the varying survival values; also, except for Phi = 1.00, individuals showed higher chances of being identified as residents than seasonal residents (Fig. 3A).

Ananias, Jesus & Yamamoto’s (2008) index was sensitive to temporality because for low values of p and Phi, it was more likely to classify individuals as temporary residents (Fig. 3B), whereas, with high probabilities of p and Phi, the probability of classifying annual resident individuals was higher. Temporary residents (i.e., occasional visitors) were surprisingly scarce even with high survival values (Phi ≈ 1.00), compared to the other definitions.

Simulations under the Martin & Silva (2004) definition of R rarely classified individuals as residents (Fig. 3C). Only at the highest survival and resighting probabilities ( ≈1.00) was it possible to identify up to 25% of the individuals as residents.

On the other hand, the definition by Möller, Allen & Harcourt (2002) was more relaxed in identifying transient and resident individuals for both low and high survival and resighting values, respectively; conversely, individuals were less prone to be classified as occasional visitors (Fig. 3D) for resighting values over 0.10, and this effect was counterintuitively higher for the maximum survival value.

Finally, simulations under the definition by Chabanne et al. (2012) had a higher chance to detect transients and occasional visitors (i.e., partial residents) (Fig. 3E) and were more conservative when classifying residents, as these occur only at p > 0.50 and Phi ≈ 1. Survival played an important role in the classification of occasional visitors since these were less likely to be identified at Phi > 0.95 and p > 0.20, turning them into resident individuals.

Validation test

The application of the R definitions to simulated data showed that in most cases different proportions of individuals in each category occurred as compared to empirical datasets (Table 3); however, patterns were similar. For instance, Martin & Silva (2004) index failed to classify individuals as “permanent residents” in both simulated and empirical data; also, Möller, Allen & Harcourt (2002) definition showed higher chances to classify resident individuals in both cases (simulated and empirical). The opposed trend was found in the case of Rosel et al. (2011) and Dinis et al. (2016). Furthermore, the indices that included occasional visitors had similar patterns for simulations and empirical data, including Ananias, Jesus & Yamamoto (2008) and Chabanne et al. (2012). Nevertheless, simulated data fitted the distribution for the Savannah dataset only under Chabanne et al. (2012) classification, as well as in Dinis et al. (2016) and Conway (2017), and these last two also fitted the Alvarado dataset (Table 3).

On the other hand, empirical results (Alvarado, Mexico and Savannah, USA) showed similar values for the three Ballance’s parameters, but the occurrence and periodicity were lower, whereas permanence was higher compared to simulations (Table 4). In this case, only the Periodicity was significantly different for both study cases.

As described earlier, from the 11 SF indices, seven were numerical, three were categorical and one had multiple categories. On the other hand, from the R indices, five were treated as categorical and three as “fuzzy” variables.

The contribution of each SF index to the total Gower distance for both empirical datasets showed a comparable pattern, where similarities among individuals were well represented in most cases (i.e., contributions >0.70) (Fig. 4). The greatest contribution corresponded to the occurrence classification by Balmer et al. (2008), the IH4 by Tschopp et al. (2018), the monthly and seasonal sighting rates by Chabanne et al. (2012), and Ballance (1990) (occurrence), whereas the least contribution corresponded to Ballance (1990) (periodicity parameter).

Conversely, the contribution of each R index to the total Gower distance for both empirical datasets were highly variable. For instance, computations for the Alvarado dataset were not possible under the Dinis et al. (2016) definition, since all individuals were classified as not resident. Moreover, contributions were well represented (i.e., contributions >0.70) for the Alvarado dataset under four of the definitions (i.e., Zolman, 2002; Ananias, Jesus & Yamamoto, 2008; Chabanne et al., 2012; Conway, 2017), but only in one instance (i.e., Chabanne et al., 2012) for the Savannah population (Fig. 5).

Performance

SF indices were not tested for performance with the simulated data, given the high redundancy in the information; that is, high co-variability was previously observed among all the indices. However, the simulations of sighting histories with the apparent survival probability of the empirical data (i.e., 0.92), showed that the R indices with the highest contributions were those of Ananias, Jesus & Yamamoto (2008), Chabanne et al. (2012), Conway (2017), and Möller, Allen & Harcourt (2002) for different resighting rates (Fig. 5), although these were highly variable. The contribution for most of the indices was within the simulated intermediate values, but the simulations showed that the classification by Chabanne et al. (2012) had more consistent results, due to less variability in the classification of the simulated resighting scenarios. Furthermore, the empirical data fell within the intermediate values of the simulated data (Fig. 5).