## Joint estimation of vital rates with a null effect, covariate effect and random effect

The previous steps have provided data on age and reproductive-status-at-death from which we can calculate vital rates. Usually, these vital rates are calculated from separate models. The Siler (1979) and Heligman and Pollard (1980) models are used to calculate vital rates associated with survival (Saavedra 2018). Logistic regression is used to calculate the vital rates associated with reproduction (Mannocci et al. 2012; Murphy et al. 2009). In the previous chapter, I described a new approach to survival estimation based on the Reed model (Reed 2011). We have developed a similar approach for estimating vital rates of reproduction by considering an Accelerated Failure Time model³ (AFT) model. Rather than considering these models separately, we combined them into a joint model.

### Cohort effect

The cohort effect is the proportion of the observed outcome that is due to the characteristics of the cohort.

Here, the cohorts are annual and have death as a common event. The cohort effect captures the effect of the common characteristics of individuals who died in the same year on the estimated vital rates.

3: Parametric model that assumes that covariate accelerate or decelerate the reaching of the event.

## Model selection

The estimation of the parameters of the two models is linked by the presence of random effects related to the year of death of the individuals. This random effect is analogous to a cohort effect (see margin-box for a definition). We considered 8 increasingly complex models. The first 2 do not involve the common cohort effect (no year effect). For each model, we computed the Widely Applicable Information Criterion (WAIC) and the Leave-One-Out cross-validation Information Criterion (LOOIC) as they allows to select the better model given the out-of-sample prediction (Vehtari, Gelman, and Gabry 2017). The better model to use is the one with the lowest WAIC or LOOIC. The selection given each model complexity is detailed in Appendix Table 1⁴. According to both WAIC and LOOIC, the better models were the 8 and the 7. From the perspective of statistical parsimony we used the model 7 as it is less complex than the 8 (no correlation between survival and reproduction estimation).

We also used the approach commonly used in the literature to estimate survival and reproduction of common dolphins (i.e. Siler, Heligman-Pollard and logistic regression). I do not detail here the functions related to these approaches, but I describe them in Appendix 7.4. Only the models and functions associated to the joint modelling approach will be detailed. Functions' parameters are estimated using a Bayesian framework. From the previously obtained data, the times to the event (death or puberty) are used to estimate the parameters of the models. These parameters are then entered into each function to calculate vital rates. The details of

4: For technical details of convergence and parameterisation see Appendix subsection 7.4

<table><caption>**Table 4.4: Models specification.**
Model 1 and 2 are not joint models. Simple: a null effect, Cov: Covariate effect, Years: Random effect (joint modelling on Years), Cor: Correlation between random effects, Trend: Trend effect on random effects. WAIC and LOOIC are depicted in Appendix Table 1.</caption><thead><tr><td>Model</td><td>Specification</td><td>Joint modelling</td></tr></thead><tbody><tr><td>1</td><td>simple</td><td>No</td></tr><tr><td>2</td><td>cov</td><td>No</td></tr><tr><td>3</td><td>years</td><td>Yes</td></tr><tr><td>4</td><td>years + cor</td><td>Yes</td></tr><tr><td>5</td><td>cov + years</td><td>Yes</td></tr><tr><td>6</td><td>cov + years + cor</td><td>Yes</td></tr><tr><td>7</td><td>cov + years + trend</td><td>Yes</td></tr><tr><td>8</td><td>cov + years + trend + cor</td><td>Yes</td></tr></tbody></table>
the estimated vital rates and the associated models are provided in Table 4.5:

<table><caption>**Table 4.5: Vital rates information.** For each important vital rate here, the mathematical notation, method of calculation and inclusion or not in the Leslie matrix model are detailed. Only age-specific survival and fecundity rate (female offspring) calculated from the joint modelling approach are used to inform the Leslie matrix. The fecundity rate is calculated using the proportion of mature females obtained from the joint modelling approach. The latter takes into account the pubescent status as a reference event.</caption><thead><tr><th>Vital rate</th><th>Notation</th><th>Calculation</th><th>Use in the Leslie matrix model</th></tr></thead><tbody><tr><td>Survivorship</td><td>l(x)</td><td>Siler & Joint model: Equation 4.2</td><td>No</td></tr><tr><td>Age-specific survival rate</td><td>S(x)</td><td>Joint model: Equation 4.3</td><td>Yes</td></tr><tr><td>Hazard rate</td><td>h(x)</td><td>Heligman-Pollard & Joint model: Equation 4.4</td><td>No</td></tr><tr><td>Proportion of mature females (with maturity status)</td><td>P<sub>M</sub>(x)</td><td>Logistic Regression</td><td>No</td></tr><tr><td>Proportion of mature females (with pubescent status)</td><td>P(x)</td><td>Equation 4.6</td><td>No</td></tr><tr><td>Fecundity rate from puberty</td><td>f<sub>P</sub>(x)</td><td>Equation 4.7</td><td>Yes</td></tr></tbody></table>

Although the approach involves joint modelling, each model is detailed separately. The join is provided by the random variable $\alpha_{it}$ (cohort effect) which is an individual effect (i) related to the year of death of the individual (t).

## Survival part

The age-at-death data obtained from the stratified random sampling are detailed in Figure 4.9. Each age within this data distribution is used to inform the parameter $x_{it}$. Age-at-death dataset from Mannocci et al. (2012) is also used to compare the effect of sampling on the survivorship value since authors did not perform stratified random sampling but used all available age-at-death values.

The survival model is the following:

$$ \log x_{it} = \mu + \sum_{j=1}^{p} \gamma_{pj} X_{ip} + \alpha_{it}^{1} + \sigma \times \epsilon_{i} - \frac{Z_{i}}{\beta} \quad (4.1) $$

It consists in modelling the logarithmic transform of $x_i$ (which is the age-at-death for each individual i) assuming a normal-Laplace distribution for $\log x_i$ where $\mu$ is a location parameter; $\sigma$ and $\frac{1}{\beta}$ are (positive) scale parameters; and $Z_i$ are independent standard exponential deviates. The parameter $\frac{1}{\beta}$ quantifies individual frailty (Kannisto 1991; Reed 2011) or persistent demographic heterogeneity (Cam, Aubry, and Authier 2016). We specified the following distribution for residuals $\epsilon_i \sim N(0, 1)$ and $\frac{1}{\beta} > 0$ which
<figure><img src="image_69.png" /><figcaption>Figure 4.9: age-at-death data as an age pyramid by sex. The maximum age class is 21-22 years with some females of this age. The data seem to be characterised by a low representation of age groups 0-1, 0-2 and sometimes 7-8 and 8-9. Over-representation is also possible for ages 9 to 14. For comparisons, see Chapter 2 with two age distributions on the same population. In total, there are 642 individuals.</figcaption></figure>

allows to accommodate a bathtub-shaped hazard rate curve which is the pattern expected for long-lived species of vertebrates in the wild (Choquet et al. 2011b). It is possible to incorporate p individual-level covariates $X_{ip}$ (such as sex or cause of death) and random variables $\alpha_t$ (such as years). The superscript $\alpha_t^1$ means that this is the random effect associated with survival model.

The vital rates associated with the survival component are:

* Survivorship $l(x)^5$ at age x with $\Phi^c$ as the cumulative density function of a standard normal distribution with $\Phi^c = 1 - \Phi$ (with $\Phi$ as the density function of a normal distribution):

$$ l(x) = \Phi^c \left( \frac{\log x - \mu}{\sigma} \right) - \exp \left( \beta(\log x - \mu) + \frac{\beta^2 \sigma^2}{2} \right) \times \Phi^c \left( \beta\sigma + \frac{\log x - \mu}{\sigma} \right) \quad (4.2) $$

* Age-specific survival $S(x)^6$:

$$ S(x) = \frac{l(x + 1)}{l(x)} \quad (4.3) $$

5: the proportion of the population that survive at least to age x.

To apply covariates and/or random and/or trend effects it is on the intercept $\mu$:

$$ \mu_{1997} = \mu + \alpha_{1997} $$

6: the proportion of the population at age x that reaches the age x + 1.
► Hazard rate $h(x)^7$ with $R$ as the Mills ratio ($R = \frac{\Phi^c}{\phi}$) which is the ratio between the cumulative density function ($\Phi^c$) and the probability density function of a standard normal distribution ($\phi$):

7: Instantaneous probability of dying at time $t+dt$ given that and individual $i$ has survived until time $t$

$$h(x) = \frac{\beta}{x} \times \frac{R \left( \beta\sigma + \frac{(\log x - \mu)}{\sigma} \right)}{R \left( \frac{(\log x - \mu)}{\sigma} \right) - R \left( \beta\sigma + \frac{(\log x - \mu)}{\sigma} \right)} \quad (4.4)$$

This approach is accurate to estimate survivorship with a sample size of 300 while for hazard rate it is better with at least 500 individuals (Rouby, Ridoux, and Authier 2021).

## Reproductive part

<figure><img src="image_70.png" /><figcaption>**Figure 4.10:** Distribution of reproductive status by age (regardless of sex). No immature individuals were founded after 11-12 years. Pubescent individuals are found from 2-3 years to 11-2 years of age with the largest proportion between 4 and 6 years. The first mature individuals are found at 4-5 years of age. For comparisons with other datasets from this population, see Chapter 2. In total, there are 240 individuals.</figcaption></figure>

The reproductive model is the following:

$$\log x_{it} = \mu + \sum_{j=1}^{p} \gamma_p X_{ip} + \alpha_t^2 + \sigma \times \epsilon_i \quad (4.5)$$

It consists in modelling the logarithmic transform of $x_{it}$ (which is age-at-death associated with a reproductive state) assuming a Weibull distribution for $x_{it}$ and a standard Gumbel distribution $\epsilon_i \sim G(0, 1)$. The Weibull distribution used is a two-parameter distribution: $\alpha$ is the shape parameter which describes the distribution pattern of the risk of becoming sexually mature; and $\lambda_i$ which is the scaling parameter defining the position of the

if $x_{it} \sim \text{Weibull}(\alpha, \lambda)$
$\log(x_{it}) \sim \text{Gumbel}(\mu, \sigma)$ with $\alpha = \frac{1}{\sigma}$
and $\lambda = e^{-\mu}$
distribution curve. Residuals $\epsilon_i \sim G(0, 1)$ which induces a Weibull law.

It is possible to include p individual-level covariates $X_{ip}$ (such as sex or cause of death) and random variables. The random variable $\alpha_t^2$ is the cohort effect associated to the reproductive model. Covariates inclusion allows to compute the proportion of mature females, which is the sex usually used in unisex animals dynamic populations, without subseting the dataset and lose information. Three reproductive status (i.e. model censoring states) are used to inform the model parameters estimation: Immature (right-censored), Pubescent (not censored, the most informative) and Mature (left-censored).

Vital rates associated with the reproductive component:

* Proportion of mature females $P(x)$ at age $x$ with $\alpha$ as the shape parameter of the Weibull distribution ($\alpha = \frac{1}{\sigma}$ with $\sigma$ as the standard deviation of the Gumbel distribution of $\log(x_{it})$) and $\mu$ as the intercept of the Gumbel distribution:

$$P(x) = 1 - \left( \exp \left( - \left( \frac{x}{\exp\left(\frac{-\mu}{\frac{1}{\sigma}}\right)} \right)^{\frac{1}{\sigma}} \right) \right) \quad (4.6)$$

<figure><img src="image_71.png" /><figcaption>Figure 4.11: Reproductive status of females given four age-class interval. The total number of female used is about 91. For a comparison with other dataset, see Chapter 2.</figcaption></figure>

* Gestation rate $G(x)$ at age $x$. The gestation rate was obtained considering 4 age-class intervals (4-8, 9-13, 14-18 and 19-22) as depicted in Figure 4.11.

Gestation rate before 4 years is equal to 0 since no female was pregnant before 4 years. It is equal to 0.125 for 4-8 since one of the 8 mature females of this age was pregnant. It is
equal to 0.45 for 9-13 since 15 of the 37 mature females of this age were pregnant. It is equal to 0.35 for 14-18 since 10 of the 29 mature females of this age were pregnant. And It is equal to 0.30 for 19-22 since 5 of the 17 mature females of this age were pregnant.

► Fertility rate $f_p(x)$ at age $x$. The proportion of mature females (estimated with the pubescent status in the AFT modelling framework) is multiplied by the gestation rate (divided by two assuming a sex ratio at birth of 1:1) and multiplied by the age-specific survival rate considering the transition to the next age $S_{x-1}$:

$$f_p(x) = P(x) \times \left( \frac{G(x)}{2} \right) \times S_{x-1} \quad (4.7)$$

The proportion of mature individuals through the ages is obtained considering both sexes. To compute fertility rate, we used the proportion of mature females computed with the pubescent status (through the covariate inclusion). It is more accurate to use females for demographic modelling since it is considered as the limiting sex in the population and it is easier to measure the reproductive output than for males.

## Leslie matrices and elasticities

The vital rates obtained above allow the construction of life tables. The tables and their construction are detailed in the supplementary material. Once vital rates have been obtained, including age-specific survival ($S_x$) and fecundity rates from pubescents ($f_{p_x}$), a Leslie matrix model can be filled in to obtain the population growth rate. This growth rate is then used in a geometric population dynamics model.

The matrix population model is the following:

$$ \begin{bmatrix} N_0 \\ N_1 \\ N_2 \\ N_3 \\ \dots \\ N_{24-1} \\ N_{24} \end{bmatrix}_{n+1} = \begin{bmatrix} N_0 \\ N_1 \\ N_2 \\ N_3 \\ \dots \\ N_{24-1} \\ N_{24} \end{bmatrix}_n \times \begin{bmatrix} f_{P_0} & f_{P_1} & f_{P_2} & f_{P_3} & \dots & f_{P_{24-1}} & f_{P_{24}} \\ S_0 & 0 & 0 & 0 & \dots & 0 & 0 \\ 0 & S_1 & 0 & 0 & \dots & 0 & 0 \\ 0 & 0 & S_2 & 0 & \dots & 0 & 0 \\ 0 & 0 & 0 & S_3 & \dots & 0 & 0 \\ 0 & 0 & 0 & 0 & \dots & 0 & 0 \\ 0 & 0 & 0 & 0 & \dots & S_{24-1} & 0 \end{bmatrix} \quad (4.8) $$

Vector $N$ gives the number of individuals in each age-class at time $t$ and then at time $t + 1$. The previously estimated vital
rates were combined into the matrix with the female age-specific survival rate $S(x)$ and the fertility rates from pubescents females $f_P(x)$. Female vital rates were obtained by including covariates without and with the random effect of year. The matrix model is therefore relative to female only. Before age class 5, the gestation rate is 0. No offspring are produced by females strictly under the age of 5.

<figure><img src="image_72.png" /><figcaption>**Figure 4.12:** Conceptual model of the Female leslie matrix modelling approach. 24 age classes are modelled. After age class 24, the females die, so there is no survival rate for this class. Before age class 5, the gestation rate is 0. No offspring are produced by females strictly under the age of 5.</figcaption></figure>

The population growth rate $\lambda$ was computed as the asymptotic growth rate of the Leslie matrix and allowed to project the remaining population at time $t$ ($N(t)$) taking into account the initial abundance at time $t_0$ ($N(t_0)$). The population's trajectory was projected through 100 years, from 2020 to 2120.

$$N(t) = N(t_0) \times \lambda^t \quad (4.9)$$

We also explored the relative changes in population's growth rate $\lambda$ caused by a proportional change in one of the life cycle parameters, referred to as elasticities (de Kroon et al. 1986). The elasticities analysis may be used to improve and study the management of long-lived species (Crouse, Crowder, and Caswell 1987).

## Baseline assumptions

* Stranded dolphin population is representative of the living population (regarding demographic features).

* Age and reproductive-status-at-death are representative of the stranded population.

* Each age-at-death is known and correctly reported.

* Each reproductive-status-at-death is known and correctly reported.
► Joint model is accurate to estimate vital rates and the survival and fecundity patterns are congruent for the species.

► The population growth geometrically and is not subject to migration.

**4.3 Results: Vital rates and population growth rate**

The first part of the results is related to the vital rates of the population. In this part, the first result highlighted is the comparison of survival estimates according to the sampling of age-at-death data. The second result highlights the patterns of survivorship and proportion of matures individuals without effect. The third result highlights the survivorship and proportion of matures individuals patterns with covariates. The fourth result highlights the survivorship and proportion of matures individuals patterns with random effects.

The second part of the result highlights the properties of the matrix model used. The first result highlights the projection of the evolution of the population (in percentage) without effect and with random effect. There is no covariate because the projection is always made from the sex covariate, to take into account the females only. The second result shows the value of the elasticities for survival and fecundity, by age and by year.

**Common dolphin population vital rates**

**Survivorship difference according to sampling**

We estimated survivorship from Mannocci et al. (2012) dataset ("All age-at-death available" Figure 4.13) and the dataset obtained in this project ("Stratified-random-sampling" Figure 4.13) using Siler (1979) and Reed (2011) models. This was done with the aim of testing the influence of stratified random sampling on the estimation of survivorship. As a reminder, Mannocci et al. (2012) did not use stratified random sampling. For the Mannocci dataset, the survivorship estimated is equal to that estimated in the Mannocci et al. (2012) study using the Siler model (Figure 4.13). Using the Reed model, there is a gap from the 0 to the 12 age-class (without taking into account the uncertainties). Then estimation are consistent and nearly overlapping. For the dataset collected in this study, both models estimated nearly the same survivorship across all age-classes even if there is a limited gap between 0 and 4 years (≈ 0.25).
<figure><img src="image_73.png" /><figcaption>**Figure 4.13:** Survivorship estimation for each age-class on the two dataset types with Reed and Siler models. The estimation realised with the Siler model is displayed in black dashed line whereas the estimation using the Reed model is displayed in black solid line. For the dataset built using all the age-at-death available (Mannocci et al. 2012), there is a difference in survivorship estimates. This difference is much less marked, or even almost non-existent, for the dataset constructed from the stratified-random-sampling.</figcaption></figure>

## Survivorship and Proportion of mature individuals with a null effect

<figure><img src="image_74.png" /><figcaption>**Figure 4.14:** Survivorship estimated from the age-at-death dataset for each age-class. The blue dotted line corresponds to the observed survivorship estimated from the Kaplan-Meier estimator. The black dashed line is made from the Reed model. Uncertainty is displayed as the 80% uncertainty interval (i.e. Bayesian credibility interval).</figcaption></figure>

The survivorship according to each age-class is highlighted in Figure 4.14. Availability of data by age is displayed in rug line. Observed survivorship from this data (obtained with Kaplan-Meier estimator) is showed in blue points. The confidence interval around the mean is a 80% confidence interval as it correspond to
the Bayesian credibility estimate. The confidence estimate is at its widest value at 10 years admitted a survivorship uncertainty of 0.3. From 0 to 10 years, the survivorship decreases from 1 to 0.25. It is almost equal to 0 at 20 years. There is a gap between the observed and the estimated survivorships from 3 to 8 years with a higher value of estimated survivorship of about 0.05.

The proportion of mature females, estimated by the joint model using pubescent status as a reference, according to each age-class is highlighted in Figure 4.16. For comparison, the proportion of mature females estimated by logistic regression from the proportions of mature females at each age is shown in Margin-Figure 4.15. In Figure 4.16, availability of data is displayed in rug line. Black rugs are associated with immature and mature individuals (both sexes). Orange rugs are associated with pubescent individuals (both sexes). Raw proportions of mature females are displayed in both Figure 4.16 (grey dots) and Margin-Figure 4.15 (black dots). The estimated fitted curve with the Accelerated Failure Time (AFT) Weibull model is expressed as a dashed line. Likewise for survivorship, the 80% interval uncertainty around the estimation is displayed in grey. The estimated ASM taking into account the pubescent biological information through the AFT approach is about 5.3 years. The maturity difference between the two approaches is therefore about 2 years.

<figure><img src="image_75.png" /><figcaption>**Figure 4.16:** Proportion of mature females estimated with joint model for each age-class. The real data are displayed as grey dots and the estimated proportion done with the AFT Weibull approach is displayed as black dashed line. The 80% uncertainty around the mean estimated is displayed as a grey area.</figcaption></figure>

<figure><img src="image_76.png" /><figcaption>**Figure 4.15:** Proportion of mature females along age-class obtained through logistic regression. Black dots are the observed proportion for each age-class. The black solid line is the model fit. Data showed a null proportion of mature females from 0 to 4 years and a full proportion of mature females from 11 to 21. From 5 to 11, the proportion increases but the pattern in the data is noisy. Female Age at Sexual Maturity (ASM) is of 7.36 years (e.g. 50% of the population is estimated as mature).</figcaption></figure>

## Survivorship and Proportion of mature individuals with Covariates

The effect of the two covariates (Bycatch and sex) on survivorship is presented in Figure 4.17 horizontal panel A. The percentage of change between the two estimated survivorship
<figure><img src="image_77.png" /><figcaption>**Figure 4.17:** Percentage change between survivorship and fertility patterns for each covariate being tested.</figcaption></figure>

A: Effect of covariates on survivorship. The difference in survival between captured and uncaptured individuals increases with age. It is about 15% at 10 years old and reaches 30% at 16 years old. Females seem to survive slightly better than males at the end of their lives, particularly from the age of 16. Their survival is about 3% better at 22 years.

B: Effect of covariates on proportion of mature individuals. The proportion of mature individuals appears to be higher in the younger age classes among the individuals caught. At 7 years of age, there are 2% more mature individuals. The females also seem to mature faster than the males. At 8 years of age, females are 5% more mature than males.

profiles is expressed in percentage. For Bycatch (or cause of death), the reference value is "Natural death". For sex, the reference value is "Male". It can be seen that the difference between profiles is more pronounced when the cause of death is taken into account in comparison to sex. By focusing on the effect of bycatch, it appears that the bycaught individuals seem to survive less well than non-bycaptured ones for all ages. The difference between bycaught and non-bycaptured individuals in survivorship goes from 0 to 10% from 0 to 8 years to 10 to 30% from 8 to 16 years.

The effect of the two covariates (Bycatch and sex) on the sexual maturity attainment is presented in Figure 4.17 horizontal panel B. The percentage of change between the two estimated sexual maturity profiles are highlighted in the same way as before. It can be seen that bycaught individuals reach sexual maturity more quickly than non-bycaptured ones. From the age of 13 to 14, however, there is no longer any difference. Sex also seems to have an effect on the estimation of the sexual maturity profile. Females seem to mature more quickly than males.
Survivorship and Proportion of mature individuals with
Random effects

<figure><img src="image_78.png" /></figure>

**Figure 4.18:** Cohort effect on survivorship. The black dashed line corresponds to the mean survivorship. Each coloured solid line of is a deviation around this mean (i.e. random effect) associated to a cohort. The cohorts share death as a common event. For example, the 2017 cohort includes all individuals who died in 2017 regardless of age. This cohort is associated with the lowest estimated survivorship. The highest estimated survivorship is for the 2002 cohort.

<figure><img src="image_79.png" /></figure>

**Figure 4.19:** Cohort effect on proportion of mature females. The black dashed line corresponds to the mean proportion of mature females. Each coloured solid line of is a deviation around this mean (i.e. random effect) associated to a cohort. The cohorts share death as a common event. For example, the 2017 cohort includes all individuals who died in 2017 regardless of age. This cohort is associated with the earliest Age of Sexual Maturity. The latest Age of Sexual Maturity is estimated for the 2002 cohort.
Mean survivorship signal (dashed black line) and the deviation from this mean signal due to the random effect of the years (coloured solid lines) are available in Figure 4.18. The year associated with the highest survivorship is 2002 whereas the year associated with the lowest survivorship profile is 2017. Years before 2007 show a profile with a higher survivorship than the mean signal. After 2007, the survivorship was associated to lower values with the lowest values recorded for recent years after 2016.
Mean proportion of mature individuals signal (dashed black line) and the deviation from this mean signal due to the random effect of the years (coloured solid lines) are available in Figure 4.19. The year associated with the latest maturity (50 % of the population estimated as mature) according to age is 2002. The year associated with the earliest maturity according to age is 2017. Years from 2009 show a profile with a latest age-at-maturity than years after 2010 that show a profile with an earliest age-at-maturity.

Population growth rate and elasticities

The mean signal of the projected population abundance
trajectory (in percent) between 2020 and 2120 from the estimated
vital rates with the covariate "female" is presented in Figure 4.20.
The mean estimated growth rate is about 0.868 which correspond
to a population's extinction within 40 to 50 years.

<figure><img src="image_80.png" /><figcaption>Figure 4.20: Remaining population from 2020 to 2120 given exponential geometric growth population model. The asymptotic growth rate is about 0.868. For a comparison with previously made population projections on this population, see Chapter 2.</figcaption></figure>

The deviations induced by the annual random effects from
the mean signal of the population abundance trajectory projection
are presented in Figure 4.21. The cohort effect (years of deaths)
deviations are displayed in colour likewise for survivorship and
proportion of mature individuals. The year with the highest growth
rate is 2002 with $\lambda = 0.928$ which corresponds to an extinction
within 70 to 80 years. 2017 is the year with the lowest growth
rate value ($\lambda = 0.848$) which corresponds to an extinction within
30-40 years. Years from 2012 show a faster extinction with much less population remaining over time compared to the years before 2012.

<figure><img src="image_81.png" /></figure>

**Figure 4.21: Cohort effect on population growth from 2020 to 2120.**
The population model admits a exponential geometric growth. The lowest asymptotic growth rate is observed for the 2017 cohort whereas the highest is observed for the 2002 cohort.

The mean proportional contribution to the growth rate, of each age-class for both survival and reproduction are displayed in Figure 4.22. For both facets, the mean signal is shown as a black dashed line and the deviations from the mean signal due to the random year effect are shown as solid coloured lines.

Regarding the survival component, the contribution decreases along the lifetime. The estimated contribution is at his highest from age 0 to 6 and then drops to 0 at age 22. From age 0 to 6, cohorts from 1997 to 2008 are associated to a lower contribution compared to the mean whereas years from 2015 to 2019 are associated with a higher contribution. From age 7 to 22 cohorts from 1997 to 2008 are associated to a higher value compared to the mean whereas cohorts from 2016 to 2019 are associated to a lower value for the same age-class interval.

Regarding the reproductive component, the contribution increases from age 1 to 8. The age-class which contributes the most to population growth rate is the $8^{th}$. Then the estimated contribution decreases up to the age-class 22. The contribution of age classes before 6 is higher than the average for the years-of-death cohorts after 2010 and lower for cohorts before 2010. Then, from age 6 until age 10 cohorts with an higher value of contribution than the mean are cohort from 2016 to 2019 with a difference peak at year 8. These years are also associated to a lower contribution compared to the mean signal from age 12 to 22. For cohorts before 2006 the contribution is lower than the mean from age 6 to 11. Then it is higher from age 12 to 22.
<figure><img src="image_82.png" /><figcaption>**Figure 4.22:** Cohort effect on elasticities of age-specific fecundity (Reproduction) and survival rates (Survival). Black dashed line corresponds to the mean signal and solid coloured lines correspond to the deviations around this mean signal. The age class whose reproduction contributes most to the growth of the population is the 8-year-old age class. The age classes whose survival contributes most to population growth are the younger age classes (0-6 years).</figcaption></figure>

Age-specific survival rates are the vital rates that contribute most to the value of the growth rate (0.9 when adding up each elasticity). age-specific fecundity rates contributes slightly (0.1 when summing each elasticity). The contribution appears to be unequal across cohorts.

The 2017 cohort has the highest contribution to population growth from reproduction and survival of young age classes (3-8 for reproduction and 0-6 for survival). It is also the cohort whose reproduction and survival of adult year-classes (14 to 22 for reproduction and 10 to 20 for survival) contribute least to population growth. This pattern is reversed for the 1997 cohort.

**4.4 Discussion: Important aspects and implication of the study**

We obtained age-at-death data for individuals (females and males) selected by stratified random sampling. We obtained reproductive-status-at-death data for all individuals (females and males) for which we had gonad samples available.

We used these data to obtain vital rates associated with survival and reproduction. The vital rates were obtained using a joint modelling approach. This approach takes into account a random effect common to both survival and reproduction, relating to the year of death of each individual.

The vital rates that will be presented in this section are presented according to the inclusion of: a null effect, a covariate effect and a random year effect.

The rates obtained by this approach were used to inform a Leslie matrix model. The values used are those for females (via the sex covariate effect) without and with random year effect. Female fecundity is obtained from modelling the attainment of puberty. The value of the gestation rate induces a female offspring equal to 0 before 5 years.
The growth rate obtained from the Leslie matrix was used to project the evolution of the population without and with random year effect. The elasticity properties of the matrix(es) were calculated.

**Representativeness bias and age distribution**

The approach is based on a cross-sectional monitoring which is fundamentally related to strong assumptions difficult to meet. The most important to meet was addressed by (Caughley 1966) which is the necessity for the population to be at a stationary age distribution when using age-at-death recorded from strandings. The stationary age distribution is an age distribution that induces a stable population whose abundance does not increase or decrease over time. This stationary age distribution is a special case of the so-called stable age distribution (Caughley 1966). The stable age distribution emerges when the birth and death rates are fixed and the distribution of ages in a population are stable (Carey and Roach 2020, p. 128). Our age-at-death dataset does not appear to be representative of a so-called stable population for a delphinid (Mannocci et al. 2012; Barlow and Hohn 1984). Since the stationary age distribution (required as a baseline assumption) is a special case of a stable age distribution, it is likely that our age-at-death sample does not meet the baseline assumption.

The age classes that could be biased are the calves (0-2 years), adult (8-12) and older-senescent (> 19 years) (Barlow and Hohn 1984; Stolen and Barlow 2003). The representation of juveniles and older-senescent individuals may be influenced by drifting processes. They may be less buoyant than adults (Peltier et al. 2012) It is therefore possible that they sink more easily. The under-representation of calves in the dataset may also be explained by a more rapid decomposition, greater vulnerability to predation, scavenging, and lower detection probability of small-sized animals (Stolen and Barlow 2003).

The Over-representation of adults in the age-at-death dataset may be explained by three assumptions (Barlow and Hohn 1984):

► 1) Non stable age distribution resulting from a disturbance.

► 2) Inconsistent deposition of GLG during life.

► 3) Oversampling resulting from a filter before stranding.

The second assumption is unlikely to be possible. For common dolphin it is commonly admitted that GLGs deposition is constantly 1 per year (Murphy and Rogan 2006). To validate or invalidate the first assumption, much more age data would be needed. It would be possible to see whether a change in representation for a given age class is repeated in subsequent years as individuals advance in age. In case of disturbance, it should propagates across years

Strandings are probably biased regarding the juvenile and the older senescent age-classes with fewer individuals than expected for a wild population.

Strandings are probably biased regarding the adults age-classes around 6 to 10 years with more individuals than expected for a wild population.
regarding age-class distribution (Barlow and Hohn 1984). The third
assumption could be possible if we consider an overmortality of
these age classes due to the bycatch. As shown in Figure 4.17, the
age classes most likely to have a drop in reduced survivorship (Age
> 8 years) when the bycatch covariate is taken into account are also
those that would appear to be over-represented in the dataset. This
observation could support the validity of the third assumption.
To better accommodate the over-representation of adults and the
under-representation of calves, it might be interesting to use a
mixture model that allows to represent two hazard rate components
(Rouby, Ridoux, and Authier 2021). The first specific to calves. The
second specific to adults. Before testing this survival model, a
simulation study should be conducted.

Contribution of stratified random sampling approach

This study is the first to employ a stratified random-sampling
framework while estimating the vital rates of a marine mammal
population from stranding data. Usually, demographic studies
from strandings use all age-at-death data available (Mannocci
et al. 2012; Barlow and Hohn 1984; Stolen and Barlow 2003). The
stratified random sampling is a way to satisfy both ecologist's and
statistician's requirements for the representative dataset (Roleček
et al. 2007). The advantage using this approach is that it allows to
maximises the probability of best representing the age structure
of the stranded population. It is important to match the stranding
signal as closely as possible on the basis of criteria (that shape
strata) influencing the demography and individual selection (such
as cause of death). Such an approach is commonly used in fisheries
studies (Holden and Raitt 1974; Cochran and William 1977; Lai
1993), forestry ecology (de Vries 1986) or human psychology and
social sciences for which it is important to correctly represent age
conditions in a population (Banerjee and Chaudhury 2010). As it is
shown in Figure 4.13 survivorship results can vary a lot depending
on the sample collected.

The sampling in the Mannocci et al. (2012) study was not
controlled and estimated survivorships may vary. The dataset
from (Mannocci et al. 2012) study were built using all the age-at-
death data available. As it was shown in Figure 4.13, the estimated
survivorship is not the same depending on the survival model
used (Siler or Reed). There is a gap between the two curves from
age 0 to about age 10. The dataset we built in this project was
constructed using stratified random sampling. The survivorship
estimates from this dataset are almost the same between the two
survival models. It is therefore possible that by not controlling for
sampling and taking all available age-at-death data, this reflects

The sampling plan seems essential to
minimise bias in estimating popula-
tion survival in cross-sectional analy-
sis.
an additional bias that may influence the estimation of survival
models. In the case where sampling is controlled, the estimated
values are the same regardless of the models used. The stratified
random sampling may help to reduce the sensitivity of the estimate
to the age distribution.

We think that it is essential to perform a stratified ran-
dom sampling (if sufficient data) in order to better represent the
stranding signal. We propose a three-steps framework to minimise
representativeness bias using strandings (also see Margin-Figure
4.4 and Diagram 4.2):

► 1) Define the need for representativeness by choosing relevant criteria and calculate the proportion of each in the stranded population

► 2) Stratify by assigning each individual a value for each criterion

► 3) Random sample individuals within each stratum to match proportions of the stranded population

This approach ensures that uncontrolled selection bias is not added
to the sample or is at least minimised. The representativeness bias
is therefore due solely to strandings.

Viability of the common dolphin population

<figure><img src="image_83.png" /></figure>

Figure 4.23: Death origin of the dolphins used in this study. The blue dots represent the age-at-death dataset. The red dots represent the reproductive-status-at-death dataset. The origin at death was obtained by using the MOTHY inverse drift model. The great majority of individuals died on the continental shelf of the Bay of Biscay.

The common dolphin population is currently considered as
a single management unit by OSPAR whose population is panmictic
(Murphy et al. 2021). In principle, the vital rates produced in this
study should refer to the whole population (Benjamins et al. 2014).
The samples used here (age-at-death and reproductive-status-at-
death) were constructed from individuals stranded on the French
coast. For some, it was possible to trace their supposed origin
at the time of their death. This is achieved through the use of
inverse drift models (Peltier et al. 2016). Almost all the individuals
died on the continental shelf of the Bay of Biscay (Figure 4.23). This does not mean that they are uniquely tied to this area. But it is a possibility to consider if we assume the hypothesis put forward by Lahaye et al. (2005) and Caurant et al. (2009) which have highlighted the possible existence of two stocks that are tied to either the continental or the oceanic part of the Bay of Biscay. If individuals were restricted to this part of European waters, then it would be difficult to generalise the vital rates obtained at the scale of the OSPAR management unit. Conversely, it would be possible to use them to inform the descriptors of the MSFD which considers parts of the population's distributional range.

It is possible that the bycatch pressure on the population is increasing and resulting in poorer survival. Since 2016, the number of stranded dolphins including bycaught individuals has been increasing almost constantly (Dars et al. 2020). The time series data we have used takes into account the year 2017. This year was marked by a very strong bycatch episode with a record number of stranded individuals (which will be exceeded by 2019 and 2020). This is the year in which survivorship is lowest, sexual maturity is reached most quickly and the population growth rate is lowest. It is therefore possible that the cohort effect observed is mainly due to an increase in bycatch pressure on this population. To compensate for this pressure, the population would tend to mature and pubesce more quickly as it is suggested by the cohort effect. The inclusion of the covariate "cause of death" may reinforce this assumption. As shown in the Figure 4.17, bycaught individuals survive less well that non bycaught ones and seem to mature more quickly. The external mortality force exerted by humans through bycatch could decrease the survival pattern of the common dolphin and speed up sexual maturation at the population level. An other interpretation of the survivorship difference between bycaught individuals and non bycaughts is sampling bias. The survival profile associated with bycaught individuals may be biased by the over-representation of bycaught individuals in the sample. The inclusion of the "cause of death" covariate would therefore limit the bias introduced by bycatch in the survival analyses.

By looking at the elasticity properties of the Leslie matrices (Figure 4.22), reproduction contributes only for 10% of the growth rate whereas survival contributes up to 90%. This differential contribution is expected for a long-lived species like the common dolphin (Heppell, Caswell, and Crowder 2000). See Margin-Figure 4.24 for a comparison with other mammals and delphinids. The low contribution of reproduction to population growth (regardless of cohort) could explain why, despite the reduction in reproductive age, the growth rate is still decreasing.

<figure><img src="image_84.png" /><figcaption>**Figure 4.24:** Stage-specific elasticities for mammal populations, grouped by age at maturity and ordered by increasing generation time.</figcaption></figure>

A. Age at first maturity = 1 yr; no juvenile stage. B. Age at first maturity = 2 yr; fertility elasticity = juvenile survival elasticity. C. Age at first maturity >2 yr. This groups correspond to species that exhibit a similar life history as the common dolphin. With delphinids in it and marine mammals. Fertility elasticities for those species are almost equals to the one that we founded for the common dolphin. Figure from Heppell, Caswell, and Crowder (2000).
Predicted but not observed extinction

Mannocci et al. (2012) projected the trajectory of the common dolphin population and predicted that the population would become extinct within 100 years at an average growth rate $\lambda = 0.945$. The average growth rate that we estimated in this study is lower with a value of $\lambda = 0.868$. Considering this rate, the population will go extinct within 50 years. This growth rate seems to vary between cohorts as highlighted by the projections taking into account random effects (Figure 4.21) but none of the years taken into account is associated with a projection as favourable as that of (Mannocci et al. 2012). Our time series covered years 1997 to 2019 whereas (Mannocci et al. 2012) used a time series from 1972 to 2006. The difference may be due to a different cohort effect between these two time series. Different characteristics or pressures experienced by our cohorts between 2007 and 2019. Considering our time series, the population's state seems worst in the recent years (2010 to 2019) with a growth rate that decreases from 0.928 in 2002 to 0.848 in 2017. According to our estimates and those of (Mannocci et al. 2012), a decrease should be observed. However, this does not seem to be the case according to the available abundance estimates (Appendix Table 1). It is possible that the population from the Bay of Biscay is a sink population receiving individuals from source stocks.

The basic assumptions made in this study may be unrealistic. Six main assumptions were made to conduct the study. The first assumption: "*Stranded dolphin population is representative of the living population*" may be unrealistic and lead to a mismatch between the growth rate obtained and the trend observed in the abundance surveys. The demographic rates of the stranded population may not be representative of the population at sea. Another baseline assumption that may account for the mismatch between the growth rate obtained and the observed abundance is the last one: "*The population growth geometrically and is not subject to migration*". This assumption suggests that the population grows or decreases geometrically (here exponentially) depending on the value of the growth rate. But the geometric population growth model used does not take into account migration phenomena which could be the cause of an apparent stationarity of the population in the Bay of Biscay.

The most favourable growth rates we have estimated are associated with the earliest years (< 2005) and are close to Mannocci et al. (2012) estimates.

The population should have been declining for years, but that is not what we are seeing. It is possible that this is a sink population. This hypothesis will be discussed in more detail in the general discussion.
# Quantifying the bycatch

5

The Common dolphin plays an important functional role in Bay of Biscay marine ecosystems as a top predator, but is threatened by bycatch, i.e. the incidental capture of non-targeted species in commercial and recreational fisheries. Bycatch has been increasing, and has been associated with a large number of winter strandings of dolphins on the French Atlantic coast since at least 2017. However, uncertainties around the true extent of common dolphin bycatch and the fisheries involved have led to delays in the implementation of mitigation measures. Current data collection on dolphin bycatch in France is with non-dedicated observers deployed on vessels for the purpose of national fisheries sampling programmes. These data cannot be considered as representative of all fisheries bycatch. This feature makes it difficult to use conventional ratio estimators, as they require a truly random sample of the fishery by dedicated observers. We applied a newly developed approach, regularised multilevel regression with post-stratification, to estimate total bycatch from unrepresentative samples and total fishing effort. The latter is required for post-stratification and the former is analysed in a Bayesian framework with multilevel regression to regularize and better predict bycatch risk. We estimated the number of dolphins bycaught for each week and ten ICES divisions from 2004 to 2020 by jointly estimating the bycatch risk, haul duration and the number of hauls per days at sea. The bycatch risk in pair trawlers flying the French flag was the highest in winter 2017 and 2019, and was associated with the longest haul durations. ICES divisions 8.a and 8.b (shelf part of the Bay of Biscay), were estimated to have the highest common dolphin bycatch. Our results are consistent with independent estimates of common dolphin bycatch from strandings. Our method shows how non-representative observer data can nevertheless be analysed to estimate fishing duration, bycatch risk and, ultimately, the number of bycaught dolphins. These weekly estimates improve current knowledge on the nature of common dolphin bycatch and can be used to inform management and policy decisions at a finer spatio-temporal scale than has been possible to date. Our results suggest that limiting haul duration, especially in winter, could serve as an effective mitigation strategy.

5.1 Introduction: About the by-catch phenomenon . . . . . 104

5.2 Methodology: Quantify the bycatch parameters . . . . . 108

5.3 Results: Estimation of variables and number of dolphins bycaught . . . . . . . . . . 119

5.4 Discussion . . . . . . . . 123
<figure><img src="image_85.png" /></figure>

Published works

This chapter is associated with three publications. The first is included as part of the introduction.
The second is not included. The third constitutes the core of the chapter:

► **1.** Gilbert, L., Rouby, E., Tew-Kaï, E., Spitz, J., Peltier, H., Quilfen, V. and Authier, M. 2021. Spatiotemporal models highlight influence of oceanographic conditions on common dolphin bycatch risk in the Bay of Biscay. Marine Ecology Progress Series. 679:195-212.
Doi: https://doi.org/10.3354/meps13894

► **2.** Authier, M., Rouby, E. and Macleod, K. 2021. Estimating cetacean bycatch from non-representative samples (I): a simulation study with regularized multilevel regression and post-stratification. Frontiers in Marine Science. 8:1459. Doi: https://doi.org/10.3389/fmars.2021.719956

► **3.** Rouby, E., Dubroca, L., Cloâtre, T., Demaneche, S., Genu, M., Macleod, K., Peltier, H., Ridoux, V. and Authier, M. Estimating bycatch from non-representative samples (II): a case study of Pair trawlers and common dolphins in the Bay of Biscay. 2022. Frontiers in Marine Science. p:1946. Doi: https://doi.org/10.3389/fmars.2021.795942.
## 5.1 Introduction: About the bycatch phenomenon

### Why does bycatch of common dolphins occur?

The structuring of multilevel trophic pelagic ecosystems is influenced by the presence of submesoscale and mesoscale oceanographic structures (Bakun 1997, 2006). These structures correspond to oceanic fronts (i.e. the frontier between two water masses of different density), upwellings (i.e. upward movement to the ocean surface of deeper cold usually nutrient-rich waters) and eddies (i.e. circular current of water). Upwellings contribute to the enrichment of the water column with nutrients through different processes, but regarding the eddies, it is the Ekman pump that is involved (Margin-Figure 5.1). In the event of a cyclonic wind (for the North Hemisphere) there will be a transfer of water from the centre to the periphery of the eddy, which will lead to an upwelling of cold water in the centre (Fieux 2019). Oceanic fronts support the enrichment of water with nutrients by the same process by allowing two different water masses to encounter and promote the retention of nutrients in a given area. If they are sufficiently stable over time, these structures allow the development of a whole trophic chain from phytoplankton, to zooplankton, to fish and finally to top-predators (Bakun 1997). Top predators, such as delphinids, integrate ecological processes at all levels of the trophic web through their dynamic distribution (Croll et al. 1998). A key element of delphinids distribution is the availability of fish, which is variable in space and time due to the oceanographic processes discussed above (Hyrenbach, Forney, and Dayton 2000).

The common dolphin is a top predator of pelagic ecosystems in the Bay of Biscay that feeds on a wide variety of fish species, some of which are targeted by fisheries (Pusineri et al. 2007; Meynier et al. 2008; Spitz et al. 2013), and whose presence depends on oceanographic conditions. The main fish species involved in dolphin-fisheries interactions are the sardine (*Sardina pilchardus*), the anchovy (*Engraulis encrasicolus*) (both are prey to the dolphin) and the seabass (*Dicentrarchus labrax*) (in trophic competition with the dolphin; Spitz et al. 2013; Figure 5.4 for an overview of the trophic chain). The latter species is targeted by the Pair trawlers (Margin-Figure 5.5), a *métiers* (i.e. the combination of gear, target species and fishing area) historically associated with bycatch of common dolphin in the Bay of Biscay (Peltier et al. 2021; as well as other fishing *métiers*). The occurrence of sea bass, and the fish on which the dolphin feeds, is indirectly conditioned by the availability of nutrients (such as phosphorus, nitrogen, potassium and iron). Nutrients have an influence on the development of phytoplankton

<figure><img src="image_86.png" /><figcaption>**Figure 5.1: Ekman pumping.** Horizontal divergence of the integrated Ekman transports gives rise to a vertical velocity at the base of the Ekman layer (i.e. "Ekman "pumping"): A. Cyclonic winds drive upwelling (Ekman suction).</figcaption></figure>

B. Anticyclonic winds drive downwelling (Ekman pumping). Figure from Colling (2001). The book Fieux (2019) also provides information on this subject.

<figure><img src="image_87.png" /><figcaption>**Figure 5.2: Drawings of fish important to the feeding ecology of the common dolphin.**</figcaption></figure>
(Bakun 1997). The presence of phytoplankton (and the diversity of species associated with it) influences the presence of zooplankton that feed on phytoplankton (Lampert et al. 1986). If the persistence and retention of nutrients are sufficient, then the development of plankton will allow predators such as fish to feed (Beaugrand and P. C. Reid 2003). Finally, large fish and dolphins feed on smaller fish (Spitz et al. 2013) and pair trawlers target fish of commercial interest (Peltier et al. 2021).

<figure><img src="image_88.png" /><figcaption>Figure 5.4: Diagram of trophic chain associated with common dolphin.</figcaption></figure>

The proximal relationship between fish, dolphins and pair trawlers could drive the spatio-temporal variability of the bycatch pattern: fish distribution is notoriously variable in both space and time (Erauskin-Extramiana et al. 2019) and is governed by many factors from dynamic oceanographic conditions. Changes in local distribution and abundance of prey species might be a substantial driver of their simultaneous presence with common dolphins, commercial predatory fish and fisheries in localized areas. In order to investigate the link between oceanographic factors and dolphin mortality at sea, we coupled biological data which are related to reverse drift trajectories obtained from the bodies of stranded dolphins with incidental capture evidence and environmental data intrinsic to various physical processes at high and low frequency obtained thanks to the modelling services of the Shom¹ (Gilbert et al. 2021). We built a mortality estimation model to study the influence of sea surface temperature, turbulence force or eddy kinetic energy (i.e. Eddies occurrence) and mean surface temperature gradient (i.e. Frontal structures occurrence and intensity) over the years 2012 to 2018. The winter and summer oceanographic processes do not appear to be linked to the same mortality processes. In winter, the Bay of Biscay environment is characterized by a seasonal cross-shore (West to East) surface temperature gradient with lowest temperature close to shore and intense frontal activity parallel to the coast (North to South) (Yelekçi et al. 2017). These frontal structures are freshwater fronts, correlated to the mixing of oceanic

1: Service Hydrographique de la Marine Nationale: https://www.shom.fr/

Two seasons of mortality with different oceanographic processes were highlighted:

Winter - Freshwater frontal structures.

Summer - Tidal frontal structures.

<figure><img src="image_89.png" /><figcaption>Figure 5.3: Scheme of pelagic pair trawler net.</figcaption></figure>

Figure from https://sng.ie/fishing/pelagic-old/.
waters and cold freshwater inputs from river plumes (Yelekçi et al. 2017). Their location and time line concur with observed patterns regarding bycatch mortality. Together with the effects of the associated oceanographic variables on the mortality index, results suggested these seasonal fronts may be targeted by both fisheries and common dolphins as areas where fish aggregate, thereby putting the latter at risk of bycatch by the former.

In July and August, the mesoscale dynamic activity of the Bay of Biscay is rather different than in winter. In summer, there are mainly fronts due to tidal flow (Yelekçi et al. 2017). Summer tidal fronts are quite consistent from one year to the next because they are correlated to a repetitive process (i.e. tides) (Yelekçi et al. 2017). In summer, the main frontal activity is a seasonal tidal front, called the Ushant front, located off western Brittany (Yelekçi et al. 2017). Its activity peaks in July and August (Yelekçi et al. 2017). The location of this typical tidal frontal structure concurred with the location of hot spots of fitted mortality index as well as with the location of bycaught common dolphin strandings, mainly occurring on the coast of the South Finistère region during the summer mortality season.

These results suggest the influence of spatio-temporal
oceanographic processes on common dolphin bycatch mortal-
ity in the Bay of Biscay. But the relationship should not be over-
interpreted as the model indirectly linked oceanographic processes
to common dolphin mortality. While this study aimed to explore
the "why" question, it is important to answer the "how much"
question. Especially in the context of the Habitat Directive and the
MSFD.

<figure><img src="image_90.png" /><figcaption>**Figure 5.5: Seasonal variation in the main hydrographic structures in the Bay of Biscay.** With: (1) winter warm slope current, (2) swaddies, (3) river plumes, (4) cold waters, (5) upwellings, (6) warm waters of Bay of Biscay, (7) slope current, (8) tidal currents. Figure from Koutsikopoulos and Le Cann (1996).</figcaption></figure>

Figure from Koutsikopoulos and Le
Cann (1996).

Bycatch observation and quantification

In France, the monitoring of cetacean bycatch in fisheries is non-dedicated (Cornou et al. 2018) (observers are non-dedicated to the bycatch observation), and the collected data are described as non representative, which makes it difficult to estimate the bycatch risk. The bycatch risk is the probability of catching at least one dolphin in a haul. From the calculation of the bycatch risk, it is possible to estimate the number of common dolphins bycaughts for a given year in a given ICES² division. This estimation is carried out mainly in the Working Group on Bycatch of Protected Species³ (WGBYC). Traditionally, using the observer data, the way to quantify the bycatch risk per Day at Sea (also refer as "bycatch rate" per day at sea) is to use a ratio-based estimator. A day at sea is defined as any continuous period of 24 hours (or part thereof) during which a vessel is present within a division and absent from port (Anonymous 2016). The number of bycaught dolphins is estimated⁴ as the

2: https://www.ices.dk/Pages/default.aspx

3: https://www.ices.dk/community/groups/Pages/WGBYC.aspx

4: The total bycatch Ŷ of species (i) by region (r) is estimated as the sum of the observed specimens (yᵢ) to observed Days at Sea (x), time total fishing Days at Sea (X) summed over ICES areas (a) of interest:

$$
\hat{Y}_{ir} = \frac{\sum y_{iar}}{\sum x_{ar}} \times \sum X_{ar}
$$
product of the ratio of the sum of observed bycaught individuals to
observed days at sea, times total fishing days at sea summed over
ICES division of interest (ICES 2019, p. 58). This ratio-estimator is a
design-based approach of the estimation of bycatch and is referred
as a Bycatch risk Assessment (BRA) by WGBYC. The accuracy of
the approach is largely dependent on the representativeness of
sampling protocol (notably the observer scheme design). Due to
the low observer coverage, the non dedicated-sampling protocol,
the non-representativeness of observer dataset and the small num-
ber of hauls with few common dolphin bycatch, estimates based
on ratios appear unreliable (Babcock, Pikitch, and Hudson 2003;
Course 2021; Authier, Rouby, and Macleod 2021). An other way to
quantify bycatch is to use strandings (Peltier et al. 2016).

The use of strandings can be referred as model-based es-
timate since it implies a modelling-framework. Here it will be
referred to as "stranding-based" estimates to better flow with the
rest of the chapter, including the discussion. This approach is
based on the drift reconstruction trajectory of the stranded indi-
vidual using the MOTHY model⁵ ("Modéle Océanique de Transport
d'HYdrocarbures") from Méteo France. Then the estimated mortal-
ity in each year is corrected by the probability of being buoyant.
This correction allows to estimate the number of dead dolphins,
included those that sank and were not found ashore (Peltier et al.
2016). Using this approach, it is not possible to discriminate the
involvement of each fishing métier in the mortality of common
dolphins. To fill this gap, the use of observer data is very important,
even if biased. Some model-based recent improvements allows to
minimise the impact of the representativeness bias. For example,
Luck et al. (2020) recently used a modelling framework to estimate
the corrected bycatch probability and then use the fishing effort
to extrapolate to the whole fishing fleet. Recently, Authier, Rouby,
and Macleod (2021) adapted the use of a the regularized multilevel
regression with post-stratification to infer more accurately bycatch
rates at finer spatio-temporal scale than it is actually done.

There is a need to obtain bycatch estimates at a fine spatial
and temporal scale and to use data from at-sea observers. In
addition to estimates of the number of dolphins bycaught, estimates
must be made for the variables associated with bycatch. The spatial
and temporal scale must be fine enough to allow for inter- and
intra-annual differences within ICES divisions. In this chapter, I
present how we quantified the bycatch risk, haul duration, number
of hauls per Days at Sea and number of bycaught dolphins from
biased observer data. We quantified these parameters for each
week from 2004 to 2020 within 10 ICES division to better represent
variations in bycatch spatio-temporality using the approach of
(Authier, Rouby, and Macleod 2021). This study involves firstly

It is impossible to use the classical
ratio-estimator approach on onboard
observer data to obtain the number
of dolphin bycaught at fine spatio-
temporal scales.

Strandings are used to obtain the
number of dolphin bycaught on a
monthly basis. The estimation is done
by inverse drift modelling. ICES esti-
mates are currently made in this way.

5: http://www.meteorologie.eu.
org/mothy/

Stranding based estimates do not al-
low the number of catches per métiers
to be estimated on a fine spatial and
temporal scale.

We used Multilevel Post-Stratification
to obtain the number of dolphins by-
caught by PTM on a weekly basis
from 2004 to 2020 in ten ICES divi-
sions. We also estimated the bycatch
risk, mean haul duration and number
of hauls per Day at Sea on the same
basis.
the analyse of historical bycatch monitoring data collected by onboard observers (from 2004 to 2020) on pair trawlers flying the French flag (hereafter called "PTM")⁶ which is called Obsmer⁷ ("*Observation des captures en mer* "). This dataset is the one considered as biased. From this dataset we estimated the bycatch risk, the haul duration and the number of hauls per Day at Sea. Then, we used a fishing effort dataset computed using the SACROIS algorithm⁸ and the median number of bycaught dolphin by haul to estimate the number of bycaught dolphins per year and ICES divisions. Both dataset were provided by IFREMER⁹ ("*Institut Français de Recherche pour l'Exploitation de la Mer* "). After presenting the involved data and the broad-lines of the applied methodology, I will compare the estimates made in this study with those of the stranding-based estimates.

## 5.2 Methodology: Quantify the bycatch parameters

In this section, the methodology used to estimate bycatch parameters (i.e. risk, duration and number of hauls) is detailed. First, a review of the materials available to perform the analyses is made. Then, a review of the methods is presented in outline.

### Materials

#### ICES divisions, fisheries and study area

The North-East Atlantic ocean can be geographically divided differently depending on the focus (for example focusing on fisheries management or on ecosystems). ICES defines some ecoregions to refer to ecosystems in its ecosystem-based management approach (EBM) and divisions to give an overview of the fisheries related aspects in these ecoregions (ICES 2021b; Figure 2.11 in Page 34 Chapter 2).

The Bay of Biscay, the English Channel and part of the Celtic seas are ecoregions associated to a high variety of fishing gears, métiers and targeted species. These areas are associated with sub-mesoscale and mesoscale oceanographic processes, such as eddies, river plumes and tidal fronts, that enhance ecosystem productivity and result in high availability of fishes, including commercial species (e.g. European seabass *Dicentrarchus labrax*, Sardine *Sardina pilchardus* or Anchovy *Engraulis encrasicolus*). Pelagic trawl fisheries, as pair-trawlers, are active within all of these areas and mainly target the anchovy, sardine, seabass and hake in divisions 8a and 8b (ICES 2020a). The study area for the estimation of common

6: Pair-Trawlers flying the French flag are referred as "PTM" for the rest of the manuscript.

7: https://sih.ifremer.fr/Ressources/Obsmer

8: https://wwz.ifremer.fr/sih_eng/Debarquements-effort-de-peche/Sacrois

9: https://wwz.ifremer.fr/
<figure><img src="image_91.png" /><figcaption>Figure 5.6: Study area in the North-East Atlantic ocean, with ICES divisions overlayed.</figcaption></figure>

dolphin bycatch by the PTM takes into account 10 ICES divisions (Figure 5.6). Depending on the bathymetry, it is possible to classify these divisions as oceanic or neritic. Divisions 7.d, 7.e, 7.f, 7.g, 7.h, 8.a, 8.b and 8.c are related to neritic ecosystems while division 7.j and 8.d are related to oceanic ecosystems.

## Datasets

Two main source of data¹⁰ were used. The first source of data that was used is a dataset called Obsmer. This dataset is collected as part of an onboard observer program set up within the Data Collection Framework of the Common Fisheries Policy and carried out by IFREMER, under the supervision of the Directorate of Fisheries and Aquaculture ("*Direction des pêches maritimes et de l'aquaculture*", DPMA). Onboarding ObsMer observers' primary duty is to register the length and weight composition of catches. Still, they have to report any bycatch if they witness such events. Obsmer data on PTM covers 4,484 hauls between 2004 and 2021, of which 82 were associated with a bycatch event of at least 1 and up to 50 common dolphins.

Obsmer provides the geographic position, timing, and du-

<footnote>10: Two dataet are used. Obsmer con-
tains the information about byctach
at the fine spatio temporal scale. It is
the biased dataset. SACROIS contains
the information about the true fishing
effort of PTMs from 2004 to 2020 on a
weekly basis. It is the population-level
data assumed as true.</footnote>
ration of hauls. Although Obsmer is aiming at a coverage of 10% and 5% of fishing effort for (level-3 métier) PTM for vessels of more than 15 meters and less than 15 meters respectively. Effort is quite low overall, ranging from 0.0 to 11% of Days at Sea (DaS; Table 5.1). The number of observed hauls with at least one bycatch record is very small: the yearly percentage of observed hauls with a bycatch event never exceeded 4.5% and was 0 in nearly half of the surveyed years. Obsmer data on PTM are an unrepresentative sample of hauls, largely because allowing an observer remains largely at the discretion of skippers (Babcock, Pikitch, and Hudson 2003; Benoît and Allard 2009).

The second dataset provides monthly estimates of total fishing effort in each division. This dataset is generated from the algorithm SACROIS developed by Ifremer and integrates data from Vessel Monitoring System, log-books and landing statistics (for boats longer than 18 meters from 1st of January 2004 and longer than 15 meters from 1st of January 2005; (Système d'Information Halieutique 2017)). The SACROIS algorithm is interesting for two main points: (1) it allows to correct errors that could exists in the integrated dataset due to recording or collecting errors ; (2) it reconstitutes métiers during fishing trip as there are not recorded in logbooks nor fish market data (Cornou et al. 2018). The SACROIS dataset provides the best available estimates of total effort, in Days at Sea (DaS), between 2004 and 2020 (Table 5.1). There are also refusals from skippers due to administrative and security reasons. Skippers must apply for permission to take observers on board and even if they decide to apply, permission may be refused for safety reasons (e.g. not enough space or rails not high enough).

These two datasets¹¹ are complementary to address the following challenge: is it possible to quantify bycatch parameters (such as bycatch risk, haul duration, number of hauls and number of dolphin bycaughts) from a non-representative sample ? The sample Obsmer contains the information on marine mammal bycatch at micro-level (i.e. on a haul resolution basis). This source of data allows to estimate the bycatch parameters (Luck et al. 2020). On the other hand, the fishing-vessels effort data from SACROIS provides the spatialized effort data at the scale of the whole fishing fleet (macro-level component). These population-level data on effort allows the post-stratification of bycatch risk estimated from observer data to obtain the number of bycaught dolphins at the level of the whole fleet (Authier, Rouby, and Macleod 2021). Descriptive statistics of both datasets are displayed in Table 5.1. Used in tandem, both datasets allow using regularized multilevel regression with post-stratification to estimate cetacean bycatch from non-representative samples (Authier, Rouby, and Macleod 2021).

<figure><img src="image_92.png" /><figcaption>Figure 5.7: Total weekly PTM fishing effort in Days at Sea for each year in division 27.8.a. The raw data suggest significant inter-annual variation but a common pattern. Total effort appears to increase between week 1 and week 28 before declining.</figcaption></figure>

<footnote>11: Obsmer provides the micro-scale information while SACROIS provides the macro-scale information.</footnote>
<caption>**Table 5.1: Descriptive statistics for Obsmer and SACROIS data displayed for each year.**</caption>

The bycatch events observed seem to be recurrent since 2007. The year associated with the highest level of observed bycatch is 2009 with 50 dolphins.

Total observer coverage of fishing effort is rarely more than 5%. The maximum PTM fishing effort covered in Days at Sea is 10.9% in 2008.

<table><thead><tr><th rowspan="2">Dataset</th><th colspan="5">Obsmer</th><th>SACROIS</th></tr><tr><th>Hauls</th><th>Average Duration (hours)</th><th>Bycatch events</th><th>Median nb of dolphins</th><th>Max. nb of dolphins</th><th>DaS (Coverage %)</th><th>Total Effort (DaS)</th></tr></thead><tbody><tr><td>2004</td><td>4</td><td>2.80</td><td>0</td><td>-</td><td>-</td><td>4 (0.0)</td><td>8 530</td></tr><tr><td>2005</td><td>5</td><td>4.26</td><td>0</td><td>-</td><td>-</td><td>4 (0.0)</td><td>8 790</td></tr><tr><td>2006</td><td>122</td><td>4.62</td><td>0</td><td>-</td><td>-</td><td>90 (1.1)</td><td>7 853</td></tr><tr><td>2007</td><td>727</td><td>3.89</td><td>6</td><td>1.5</td><td>5</td><td>401 (6.4)</td><td>6 305</td></tr><tr><td>2008</td><td>554</td><td>4.81</td><td>6</td><td>1.5</td><td>4</td><td>328 (10.9)</td><td>3 011</td></tr><tr><td>2009</td><td>464</td><td>5.50</td><td>20</td><td>2</td><td>50</td><td>326 (7.4)</td><td>4 413</td></tr><tr><td>2010</td><td>305</td><td>3.52</td><td>1</td><td>4</td><td>4</td><td>159 (3.5)</td><td>4 486</td></tr><tr><td>2011</td><td>173</td><td>3.99</td><td>2</td><td>3</td><td>3</td><td>86 (2.1)</td><td>4 001</td></tr><tr><td>2012</td><td>210</td><td>3.58</td><td>4</td><td>4</td><td>8</td><td>96 (2.4)</td><td>4 005</td></tr><tr><td>2013</td><td>128</td><td>3.81</td><td>2</td><td>5.5</td><td>9</td><td>75 (1.8)</td><td>4 192</td></tr><tr><td>2014</td><td>114</td><td>4.44</td><td>0</td><td>-</td><td>-</td><td>78 (1.9)</td><td>4 136</td></tr><tr><td>2015</td><td>136</td><td>2.77</td><td>1</td><td>2</td><td>2</td><td>78 (1.7)</td><td>4 597</td></tr><tr><td>2016</td><td>156</td><td>4.75</td><td>5</td><td>3</td><td>10</td><td>106 (2.3)</td><td>4 603</td></tr><tr><td>2017</td><td>196</td><td>5.23</td><td>12</td><td>2</td><td>20</td><td>124 (2.6)</td><td>4 835</td></tr><tr><td>2018</td><td>184</td><td>3.85</td><td>1</td><td>1</td><td>1</td><td>102 (2.8)</td><td>3 613</td></tr><tr><td>2019</td><td>438</td><td>5.45</td><td>11</td><td>2</td><td>8</td><td>289 (7.4)</td><td>3 139</td></tr><tr><td>2020</td><td>123</td><td>3.69</td><td>2</td><td>2</td><td>3</td><td>70 (4.0)</td><td>1 686</td></tr></tbody></table>

## Methods

The modelling framework involved two steps. First, it is necessary to both estimate the bycatch risk and the duration of hauls for each haul *i*, ICES division *j*, week *t* and year *k*. Also, it is possible to estimate the number of hauls per day at sea at the week, division and year scale. This is done with the Obsmer dataset. Then, it is necessary to rescale these estimations at the entire fleet level to estimate the number of bycatch events and number of bycaught dolphins for the same levels. This post-stratification is done with the SACROIS dataset.

Estimations were carried out in a Bayesian framework using programming language Stan (Carpenter et al. 2017b) called from R v.4.0.1 (R Core Team 2020) with library Rstan (Stan Development Team 2020). Stan uses Hamiltonian dynamics in Markov chain Monte Carlo (MCMC) to sample values from the joint posterior distribution (Carpenter et al. 2017b). Four chains were initialised from diffuse random starting points and run for a total of 2,000 iterations, discarding the first 1,000 as warm-up. Default settings for the No-U-Turn Sampler (NUTS) were changed to 0.99 for adapt delta and 15 for max tree depth (Hoffman and Gelman 2014). We fitted a total of 6 models of differing complexity (see Margin-Table 5.2 for specification and article in Appendix 7.4 for model selection): we compared models assuming either a gamma or a log-normal likelihood for haul duration, and models assuming additive effects

<table><caption>**Table 5.2: Models specification.**</caption><thead><tr><th>Model</th><th>likelihood for duration</th><th>Specification</th></tr></thead><tbody><tr><td>•M<sub>6</sub></td><td>gamma</td><td>ICES division × week × year</td></tr><tr><td>•M<sub>5</sub></td><td>log-normal</td><td>ICES division × week × year</td></tr><tr><td>•M<sub>4</sub></td><td>gamma</td><td>ICES division + week × year</td></tr><tr><td>•M<sub>3</sub></td><td>log-normal</td><td>ICES division + week × year</td></tr><tr><td>•M<sub>2</sub></td><td>gamma</td><td>ICES division + week + year</td></tr><tr><td>•M<sub>1</sub></td><td>log-normal</td><td>ICES division + week + year</td></tr></tbody></table>
versus interactive effects of week, year and divisions. Model fitting was carried out on the supercomputer facilities of the "Mésocentre de calcul de Poitou Charentes (Université de Poitiers/ISAE-ENSMA/La Rochelle Université)". Codes are available at https://gitlab.univ-lr.fr/mauthier/cdptmbycatch.

## Modelling Bycatch risk, duration of hauls and number of hauls

The objective of this part of the modelling framework is to model the bycatch risk and the haul duration at the haul level and the number of hauls per day at sea using only the Obsmer dataset. The problem with trying to estimate the value of these parameters is that we do not know how to interpolate between the observed points. Lets start with an example. In the Margin-Figure 5.8, we see that we only have 6 data points (for weeks *t* = 3, 4, 10, 24, 46 and 53). This number is reduced to 4 for a given year *k* and is reduced to 3 for a given division *j* in a given year *k* (Figure 5.10 to see the data depletion associated with stratification). We do not know how to fit a line between those points. And even less is known about how to extend the estimation of these profiles beyond the weekly scale to estimate by week for each year in each division. We only have some observations at the haul level scale *i* for the bycatch risk and the haul duration and we also have the number of hauls. Obsmer data are sparse (Margin-Table 5.3) and the stratification of these data increases sparseness (Figure 5.9).

<figure><img src="image_93.png" /><figcaption>**Figure 5.8: Diagram of bycatch risk per weeks.** This diagram shows how sparse the data available on the collection by observers is.</figcaption></figure>

<figure><img src="image_94.png" /><figcaption>**Figure 5.9: Diagram of stratification.** This diagram shows the stratification of the data. The dots represent trawl hauls. Red dots are those with at least one bycatch event. The blue dots are those with no bycatch events. The lines are stratified by week, division and year.</figcaption></figure>

To interpolate between the observed data points, we used firstly a random walk and then Gaussian Process. See Görtler, Kehlbeck, and Deussen (2019) for an interactive web-page explanation of Gaussian Processes. The vector **e** aggregates the mean weekly effects (on the linear predictor scale, e.g. logit scale for
<table><caption>Table 5.3: Dummy dataset to illustrate the typical sparseness in Obsermer data.</caption><thead><tr><th>Haul (i)</th><th>Duration (hours)</th><th>Bycatch</th><th>Week (t)</th><th>Year (k)</th><th>Division (j)</th></tr></thead><tbody><tr><td>1</td><td>1.5</td><td>No</td><td>2</td><td>1</td><td>2</td></tr><tr><td>2</td><td>5</td><td>Yes</td><td>2</td><td>1</td><td>2</td></tr><tr><td>3</td><td>5</td><td>No</td><td>5</td><td>10</td><td>8</td></tr></tbody></table>

bycatch risk) which are modelled with a first order random walk to ensure some smoothness in between-week variations (Figure 5.10; Authier, Rouby, and Macleod 2021).

<figure><img src="image_95.png" /><figcaption>**Figure 5.10:** Diagram of the estimation process. This diagram is based on the example of the bycatch risk. It could also be based on the other two estimated parameters, which are: The Haul and the number of Hauls per Days at Sea. The raw bycatch risk data are logistically transformed (logit function) in order to place the estimate on a continuous scale and not between 0 and 1. Each data point corresponds to a haul. A first estimate of the risk in one year (week 1 to 53) is modelled by a random walk of order 1. The first stratification is performed and the value of each risk profile for each year *k* is estimated from a Gaussian process whose mean is the previous weekly risk profile. The second estimation then takes place and the value of the risk profile for each division *j* within each year *k* is estimated by another Gaussian process whose mean is the previous Gaussian profile.</figcaption></figure>

From the observed data points, it is possible to infer the bycatch risk profile for each strata using a Bernoulli model (i.e. logit scale; Figure 5.10). As it is depicted below on Equations 5.1, the intercept of the fitted curve corresponds to $\epsilon_t$ and is equal to $\mu$.
The first order random walk then allows to estimate the following weekly value with a normal distribution: $\epsilon_{t+1} \sim N(\epsilon_t, \sigma_{\text{week}})$ and so on until the 53rd week.

From this overall profile (i.e. the vector $\epsilon$), a first Gaussian Process allows to estimate the bycatch risk profile for year $k$, $\delta_k$. Conceptually, the idea is to model $\delta_k$ as a deviation from the mean profile $\epsilon$. This can be described as modelling a random effect except that curves rather than scalars are being modelled here. The profile is estimated with the Gaussian Process as follows: $\delta_k \sim \mathcal{GP}(\epsilon, \Sigma_\delta)$ with $\Sigma_\delta$ as the covariance matrix.

From this year-specific profile (i.e. the vector $\delta_k$), a second Gaussian Process allows to estimate the profile of $\alpha_{jk}$ which depends on the year $k$ and on the division $j$ (Figure 5.10). The profile is estimated with the Gaussian Process as follows: $\alpha_{jk} \sim \mathcal{GP}(\delta_k, \Sigma_\alpha)$ with $\Sigma_\alpha$ as the covariance matrix.

$$ \left\{ \begin{array}{ll} \alpha_{jk} & \sim \mathcal{GP}(\delta_k, \Sigma_\alpha) \\ \delta_k & \sim \mathcal{GP}(\epsilon, \Sigma_\delta) \\ \epsilon_t & = \mu \quad t = 1 \\ \epsilon_{t+1} & \sim N(\epsilon_t, \sigma_{\text{week}}) \quad t > 1 \end{array} \right. \qquad (5.1) $$

Simpler models without such interactions, and with only additive effects, were also fitted to the data. The simplest model included only additive random (unstructured) effects (dropping the superscript for convenience):

$$ \left\{ \begin{array}{ll} \alpha_{jkt} = \epsilon_t + \delta_k^* + \alpha_j^* & \\ \alpha_j^* \sim N(0, \sigma_{\text{division}}) & \forall j \\ \delta_k^* \sim N(0, \sigma_{\text{year}}) & \forall k \\ \epsilon_t = \mu & t = 1 \\ \epsilon_{t+1} \sim N(\epsilon_t, \sigma_{\text{week}}) & t > 1 \end{array} \right. \qquad (5.2) $$

To simplify the equations, I did not include the superscripts in the Equations 5.1 but there are three alphas $\alpha^1, \alpha^2, \alpha^3$ and three deltas $\delta^1, \delta^2, \delta^3$. The superscript 1 refers to the bycatch risk, the 2 to the haul duration and the 3 to the number of hauls.

The Bycatch risk $p_{jkt}$ is modelled as follow with $y_{ijkt}$ corresponding to bycatch event (0 or 1) of the $i^{th}$ haul for each ICES statistical division $j$ in week $t$ of year $k$. The occurrence of a bycatch event is modelled as the logit transform of $\alpha_{jkt}^1$ in a binary outcome, provided by the Bernoulli law. Depending on the bycatch risk, a

<figure><img src="image_96.png" /><figcaption>Figure 5.11: Diagram of the interaction between the equations in the estimation. The value of $\mu$ gives the intercept for the weekly risk profile. The vector $\epsilon$ corresponding to the weekly risk profile.</figcaption></figure>

This vector is the mean estimate of the first Gaussian process for each year $k$. The estimate provides a risk profile for each year $\delta_k$.

From this annual vector, a final Gaussian process takes it as an average to estimate the profile of each division $j$ within each year $k: \alpha_{jk}$.

The values of this vector are then used to inform the estimation of the bycatch parameters.
bycatch event can occur as:

$$ 
\begin{aligned}
y_{ijkt} & \sim \text{Bernoulli}(p_{jkt}) \\
& \sim \text{Bernoulli}(\text{logit}^{-1}(\alpha_{jkt}^{1}))
\end{aligned} 
\qquad (5.3) 
$$

The fishing duration $d_{ijkt}$ at the haul level $i$ is modelled for the same levels $jkt$ with a Gamma distribution $\mathcal{G}$ of shape parameter $\beta$ and scale parameter $\frac{\beta}{d_{jkt}}$. The parameter $\bar{d}_{jkt}$ is associated to the second alpha as: $\bar{d}_{jkt} = e^{\alpha_{jkt}^2}$. This law allows the representation of continuous non negative value which is appropriate for a fishing duration.:

$$ 
\begin{aligned}
d_{ijkt} & \sim \mathcal{G}\left(\beta, \frac{\beta}{\bar{d}_{jkt}}\right) \\
& \sim \mathcal{G}\left(\beta, \frac{\beta}{e^{\alpha_{jkt}^2}}\right)
\end{aligned} 
\qquad (5.4) 
$$

and the number of hauls $n_{jkt}$ for the same levels $jkt$ per Days at Sea $DaS_{jkt}$ is modelled assuming a zero-truncated Poisson likelihood (by definition, a DaS implies at least one haul). The parameter $\lambda_{jkt}$ is associated to the third alpha as: $\lambda_{jkt} = e^{\alpha_{jkt}^3}$. This distribution allows the representation of non negative and non null discrete events number which is appropriate for a number of fishing hauls.:

$$ 
\begin{aligned}
n_{jkt} & \sim \mathcal{P}^{+}(\text{DaS}_{jkt} \times \lambda_{jkt}) \\
& \sim \mathcal{P}^{+}(\text{DaS}_{jkt} \times e^{\alpha_{jkt}^3})
\end{aligned} 
\qquad (5.5) 
$$

The smoothness in $\alpha_{jk}$ and $\delta_k$ is controlled via the covariance matrices $\Sigma_\alpha$ and $\Sigma_\delta$ respectively as:

$$ 
\left\{ 
\begin{array}{l}
\Sigma_{\alpha} = \Delta_{\alpha} \Omega \Delta_{\alpha}, \quad \Delta_{\alpha} = 
\begin{bmatrix}
\sigma_{\text{zone}} & 0 & 0 \\
0 & \ddots & 0 \\
0 & 0 & \sigma_{\text{zone}}
\end{bmatrix} \\
\\
\Sigma_{\delta} = \Delta_{\delta} \Omega \Delta_{\delta}, \quad \Delta_{\delta} = 
\begin{bmatrix}
\sigma_{\text{year}} & 0 & 0 \\
0 & \ddots & 0 \\
0 & 0 & \sigma_{\text{year}}
\end{bmatrix}
\end{array}
\right. 
\qquad (5.6) 
$$

Equations 5.1 and 5.6 allow to model an interaction between week, year and division. The joint model defined in Equations 5.3,

<figure><img src="image_97.png" /><figcaption>Figure 5.12: Diagram of the deviation when estimating the profile of the parameter of interest. What causes the profile to deviate is the variance due to the year effect ($\sigma_{\text{year}}$) and the variance due to the division effect ($\sigma_{\text{zone}}$). These two variances are included in the variance covariance matrices $\Delta$.</figcaption></figure>
5.4 and 5.5 includes a time-varying component at the week-scale with an interaction with year and division.

As it is showed in Equation 5.6, there are two important component to build the covariance matrix. The first is $\Delta$ which corresponds to the variance matrix either for $\alpha$ (diagonal of $\sigma_{\text{zone}}$) or for $\delta$ (diagonal of $\sigma_{\text{year}}$). These $\sigma$ are controlling the magnitude of variations from a mean profile (Margin-Figure 5.12 for a diagram).

The second is $\Omega$, which is the correlation function between the variance matrixes. It is the correlation between week $t$ and week $t'$. This function is a Matérn correlation function of order $\frac{3}{2}$ and range parameter fixed to $\frac{3}{2}$ as: $\Omega(t, t') = \left(1 + \frac{2\sqrt{3}\times d(t-t')}{3}\right) \times \exp\left(-\frac{2\sqrt{3}\times d(t-t')}{2}\right)$ where $d(t-t') = |t-t'|$ is the temporal distance (in weeks) between weeks $t$ and $t'$. The choice of the range parameter induces a temporal correlation of 0.05 after 4 weeks (Margin-Figure 5.13 and Figure 5.14; Authier, Rouby, and Macleod 2021).

<figure><img src="image_98.png" /><figcaption>**Figure 5.13:** Diagram of autocorrelation independance. The omega-covariance function admits independence after 4 weeks.</figcaption></figure>

<figure><img src="image_99.png" /><figcaption>**Figure 5.14:** Graphical representation of the assumed correlation matrix $\Omega$. The correlation is equal to 0 after 4 weeks.</figcaption></figure>

The methodology is shown in the Directed Acyclic Graph Figure 5.15 taking for example the bycatch risk estimation. We have information about bycatch risk from the dataset Obsmer. From the estimations of hyperparameters (notably $\mu$ and $\sigma_{\text{week}}$) we can estimate the others. Finally it is possible to fit the bycatch risk profile (such as the haul duration and the number of hauls per Day at Sea) with a Bernoulli model (as in equation 5.3).
<figure><img src="image_100.png" /><figcaption>**Figure 5.15:** Directed Acyclic Graph of the methodology. The observed data are the basis of the statistical inference. Each parameter can be a hyperparameter of the others. For example, μ is a hyperparameter of the whole process while α is a simple parameter. Each parameter, or hyperparameter, is estimated from the previous one and a small deviation σ induced by stratification.</figcaption></figure>

Modelling assumptions

► 1. Bycatch risk, fishing duration and number of hauls per Days at Sea may be correlated. The correlation is induced at the week level via parameter α.

► 2. There is a temporal Independence after 4 weeks.

Estimating the total number of hauls and bycatch events

The goal of this part of the modelling framework is to scale-
up the estimated bycatch risk and the haul duration at the haul
level and the number of hauls per day at sea at the fleet level using
the fishing effort dataset (Days at Sea). It allows to estimate the
total number of unobserved hauls N̂$_{jkt}$ and the number of bycatch
events Bycatch$_{jkt}$ for each week t, year k and division j.

The number of unobserved hauls N<sub>jkt</sub> that happened in
ICES statistical division *j* in week *t* of year *k* can be estimated from
the number of observed DaS in Obsmer (DaS<sub>jkt</sub><sup>Obsmer</sup>) and from total
effort DaS<sub>jkt</sub><sup>tot</sup> (and accounting for zero-truncation):

$$
\hat{N}_{jkt} = \frac{\lambda_{jkt}}{1 - e^{-\lambda_{jkt}}} \times (\hat{\text{D}}\text{aS}_{jkt}^{\text{tot}} - \text{DaS}_{jkt}^{\text{Obsmer}}) \quad (5.7)
$$
The total number of bycatch events in ICES statistical division *j* in week *t* of year *k* is estimated as the sum of events observed in Obsmer (Bycatch<sub>jkt</sub><sup>Obsmer</sup>) and the number of unobserved hauls multiplied by bycatch risk (*p*<sub>jkt</sub><sup>hat</sup>):

$$ \text{Bycatch}_{jkt} = \text{Bycatch}_{jkt}^{\text{Obsmer}} + \hat{N}_{jkt} \times \hat{p}_{jkt} \quad (5.8) $$

Similarly, for each year, the number of common dolphins bycaught in PTM can be estimated using the observed number of bycaught dolphins in Obsmer, the estimated number of unobserved hauls (Eq. 5.7), bycatch risk and either the median number of dolphins involved in a bycatch event (Table 5.1, or the grand median of *m* = 2 for years with no observed bycatch event). We used the median to attenuate the influence of some bycatch events involving up to 50 dolphins (Table 5.1). These estimates are thereafter referred to as model-based estimates.

## Baseline assumption for post-stratification

The total fishing effort data are accurate.

## Comparison with strandings

One difficulty while estimating the number of bycaught dolphins is that it is difficult to ground-truth. The sample provided by Obsmer, a non-dedicated observer scheme of marine mammal bycatch, may not be representative of all bycatch. In addition it provides very sparse data, with less than 100 observed events over 17 years when strandings have reached several hundreds per week in recent years (ICES 2020b) (for all causes of death). The number of stranded common dolphins with evidence of bycatch can be used to estimate the total bycatch mortality with reverse drift modelling (Peltier et al. 2016) used in international working groups (ICES 2020b). Reverse drift modelling corrects for at-sea drift conditions, but cannot inform on which fishing gears was responsible for bycatch¹². Hence, strandings-based estimates are considered here as total estimates of bycatch, and can be used to partially ground-truth model-based estimates of bycatch by PTM.

These model-based estimates use data independent from strandings, but they should not exceed stranding-based estimates. Secondly, whether model-estimates correlate with stranding-based ones is of interest to shed light on the increased mortality reported in the Bay of Biscay (Peltier et al. 2021). For each year we checked the magnitude of model-based estimates against stranding-based ones, and computed Pearson's correlation coefficient between the

<footnote>12: While it is possible to spot bycatch evidence as it is shown in Figure 2.19 in Chapter 2, it is very difficult to identify fishing gear involved. Sometimes it is possible to see the mesh marks of the net and assume the type of net. But this remains very exploratory and qualitative.</footnote>
two time-series at the month level. To account for drift, these
correlation were computed with and without a lag of 2 weeks
when aggregating model-based estimates at the month level.

5.3 Results: Estimation of variables and number of dolphins bycaught

Bycatch risk, haul duration and number of hauls per
Days at Sea

<figure><img src="image_101.png" /><figcaption>**Figure 5.16: Model-based estimates (posterior medians) regardless of division.** Estimates were made for bycatch risk, haul duration and number of hauls per DaS of pair-trawlers flying the French flag operating in the study area. Each colour represents a different year (parameters $\alpha_{jkt}$ in Eq. 5.1) and the dotted black line the yearly average (parameters $\epsilon_t$ in Eq. 5.1).</figcaption></figure>

Haul duration, hauls per DaS and bycatch risk per haul (Eqs. 5.7 and 5.8) were jointly estimated. Their temporal variations are displayed on Figure 5.16 for each week between 2004 and 2020. Haul duration was the highest in week 1 with a posterior median estimate of 5.8 hours that decreased to 4.0 h in week 16, before dropping to 2 h in week 24. Haul duration increased up to 3 hours in week 32 and plateaued until the end of the year. Remarkable years were 2008 and 2014 with the lowest haul durations estimated
from week 1 to 24. Other remarkable years were 2017, 2019, and 2020 with the longest haul durations estimated from week 1 to 10. In 2016, an increase in haul duration was estimated in week 48 (5 h vs 3 h on average across years).

Bycatch risk was maximum in week 1 (around 0.1) and decreased to almost 0 from week 8 onwards. 2017, 2019 and 2020 were the years with the highest estimated bycatch risk in the first 8 weeks. In particular, risk was as high as 0.20 in 2017 for the first four consecutive weeks of the year. Two years prior to 2012 were associated with an increased risk between weeks 30 and 36. Year 2016 showed a rise in bycatch risk in week 48. The weeks in the highest risk years were also those with the longest hauls. The estimated number of hauls per DaS varied little from year to year, but varied substantially within a year.

<figure><img src="image_102.png" /></figure>

**Figure 5.17: Temporal patterns for division 27.8.a.** The temporal pattern is displayed for 4 variables: haul duration (hours), bycatch events (n), bycatch risk per haul (probability), common dolphins bycaught (n). The estimated haul duration is the highest from week 1 to 20. The lowest at week 24 and then it remains stable from 32 to the end.

The estimated bycatch risk is highest between weeks 1 and 8 and is highest for years beyond 2016.

<figure><img src="image_103.png" /></figure>

**Figure 5.18: Temporal patterns for division 27.8.b.** The temporal pattern is displayed for 4 variables: haul duration (hours), bycatch events (n), bycatch risk per haul (probability), common dolphins bycaught (n). The estimated haul duration is the highest from week 1 to 20. The lowest at week 24 and then it remains stable from 32 to the end.

The estimated bycatch risk is highest between weeks 1 and 8 and is highest for years beyond 2016 and the year 2012.
There were noteworthy differences between divisions regarding bycatch risk (Figures 5.17 and 5.18). The overall signal was similar to the one observed on Figure 5.16 with the highest risk values estimated between weeks 1 and 8. Risk in 2017 and 2019 was higher by a factor 5 in week 1 compared to other years. After week 8, this difference disappeared. With respect to divisions, division 8.a. was the one with the highest bycatch risk, with an estimate as high as 0.50 in winter 2017 and 2019. For division 27.8.a, estimated bycatch events were more numerous from week 1 to 8 with high peaks for years > 2016. There was also a peak for year 2016 in week 48. The estimated number of common dolphins presented the same pattern as the number of bycatch events. For division 27.8.b, estimated bycatch events were more numerous from week 1 to 12 with a high peak for 2009 and 2012. There were two peaks for 2016, one in week 4 and one in week 10. The estimated number of common dolphins presented the same pattern as the number of bycatch events.

## Number of Bycaught dolphins

<table><caption>**Table 5.4:** Model-based estimates of common dolphin bycaught in PTM in the study area. Divisions 7.j and 8.d are labelled 'Oceanic', divisions 7.defgh are labelled 'Neritic 7' and divisions 8.abc are labelled 'Neritic 8'. Estimates (posterior median) are reported with the lower and upper bound of a 80% credibility interval (Louis and Zeger 2009)</caption><thead><tr><th>Year</th><th>Neritic 7</th><th>Neritic 8</th><th>Oceanic</th><th>Total</th></tr></thead><tbody><tr><td>2004</td><td>0 48 248</td><td>0 177 876</td><td>0 1 10</td><td>0 227 1134</td></tr><tr><td>2005</td><td>0 56 302</td><td>0 235 1101</td><td>0 2 15</td><td>0 293 1417</td></tr><tr><td>2006</td><td>0 77 378</td><td>0 208 923</td><td>0 0 3</td><td>0 286 1303</td></tr><tr><td>2007</td><td>15 45 102</td><td>0 29 111</td><td>1 2 8</td><td>16 77 219</td></tr><tr><td>2008</td><td>1 18 63</td><td>11 46 125</td><td>0 0 4</td><td>12 65 190</td></tr><tr><td>2009</td><td>10 94 248</td><td>172 315 568</td><td>0 1 6</td><td>183 412 820</td></tr><tr><td>2010</td><td>0 119 537</td><td>4 112 454</td><td>0 0 3</td><td>4 232 994</td></tr><tr><td>2011</td><td>9 128 359</td><td>0 61 270</td><td>0 1 7</td><td>9 191 635</td></tr><tr><td>2012</td><td>22 233 667</td><td>0 129 511</td><td>0 3 13</td><td>23 366 1190</td></tr><tr><td>2013</td><td>13 315 1086</td><td>0 105 442</td><td>0 5 25</td><td>13 426 1552</td></tr><tr><td>2014</td><td>0 33 158</td><td>0 50 224</td><td>0 0 3</td><td>0 84 384</td></tr><tr><td>2015</td><td>0 14 71</td><td>2 78 368</td><td>0 1 8</td><td>2 94 446</td></tr><tr><td>2016</td><td>0 15 76</td><td>55 255 852</td><td>0 0 3</td><td>55 270 929</td></tr><tr><td>2017</td><td>0 18 61</td><td>156 600 1355</td><td>0 0 1</td><td>156 618 1415</td></tr><tr><td>2018</td><td>0 2 15</td><td>1 31 147</td><td>0 0 2</td><td>1 35 163</td></tr><tr><td>2019</td><td>0 12 40</td><td>59 203 391</td><td>0 1 6</td><td>59 216 441</td></tr><tr><td>2020</td><td>0 6 27</td><td>4 50 159</td><td>0 0 5</td><td>4 57 190</td></tr></tbody></table>

The estimated total number of bycaught dolphins for each year is reported in Table 5.4. The study area was further divided into three strata: a neritic stratum in ICES subarea 7 (divisions 7.defgh) and another in subarea 8 (division 8.abc); and an oceanic stratum spanning subareas 7 and 8 (divisions 7.j and 8.d). Estimates were the lowest in the oceanic stratum of the study area, and the largest in the neritic stratum spanning ICES subarea 8. The largest bycatch estimate was in 2017, with a posterior median of > 600 common dolphins bycaught in PTM operating in the neritic stratum
spanning ICES subarea 8. There were large between-year variations
in estimates, ranging from less than a hundred (in 2018) to more
than one thousand (in 2017). Uncertainties around model-based
estimates were also large. Weekly profiles for the estimated number
of common dolphins is also available for both divisions 8.a and
8.b in Figures 5.17 and 5.18.

Comparison and correlations with strandings

<figure><img src="image_104.png" /><figcaption>**Figure 5.19:** Comparison of model-based (orange solid line) and stranding-based estimates (dark violet solid line) of **common dolphin bycatch**. Upper panel: raw estimates with uncertainty intervals (80% for model-based estimates and 95% for stranding-based estimates); middle panel: standardized (mean centered and unit variance) estimates. Bottom panel: year-level and month-level (within each year) correlation were computed. At the within year (between month level), correlations were computed with and without a two-weeks lag. Change in the magnitude of the Pearson correlation coefficient is shown with an arrow.</figcaption></figure>

Strandings data were used to estimate common dolphins
mortality due to fisheries (Peltier et al. 2016) for each month from
1990 to 2020. Stranding-based estimates aggregate mortality due
to all fisheries and do not distinguish between gears or métiers.
Nevertheless, we correlated stranding-based estimates with our
model-based estimates of mortality from PTM flying the French
flag within each year (Figure 5.19). Correlations were computed on
raw and standardized (mean centered and unit variance) values
(Figure 5.19).

Model-based estimates of bycatch by PTM were always be-
low stranding-based estimates (which do not allow to disaggregate
by *métiers*) except for 2010 (Figure 5.19). In 2010, model-based and
stranding based estimates were 465 and 343 with a large overlap
in credibility interval. Correlations between the two time-series
were always positive. The estimated correlation decreases from
0.45 in 2004 to 0.35 in 2020 (profile without lag). The estimated
correlation decreases from 0.75 in 2004 to 0.5 in 2020 (profile with
lag; Figure 5.19). The temporal trend in within year correlation was
negative over the study period.

**5.4 Discussion**

Variations of bycatch risk and number of bycaught
dolphins

The estimated bycatch risk pattern shows some variations be-
tween and within years. This risk directly influences the estimated
number of bycaught dolphins which also displays some variations.
Here, a discussion about these variations for both components is
proposed.

Inter-Year variations

The variation of bycatch risk, and subsequently of the num-
ber of bycaught dolphins, between years may be linked to ecological
factors. Bycatch can be viewed as the resultant of two components
that are the fishing activity and the dolphin presence (Margin-
Figure 5.20). The oceanographic processes may influence the two
components to enhance the encounter probability between fish-
ing engines such as PTM and common dolphins (Gilbert et al.
2021). The highest estimated bycatch risk values were reported
for the eight or ten first weeks of each year within each division
of the study area 5.16. We can speculate that the years associated
with a high bycatch risk were also those of notable spatial overlap
favoured by particular oceanographic processes. Fishing activities
and their operational conditions can also influence bycatch risk.
For example, the setting of quotas on particular fish may cause
operational changes in certain seasons. Market effects can also have
an impact on fishing activity. If a particular fish is too devalued,
it will no longer be targeted. These aspects will not be discussed
in detail here as they are multifactorial and are more a matter of
fisheries management.

Variations in the number of bycaught dolphins between years can

<figure><img src="image_105.png" /><figcaption>Figure 5.20: Simple diagram of by-catch phenomenon components.</figcaption></figure>
also be explained by other factors than the bycatch risk. They are model-based estimates obtained through the bycatch risk, total fishing effort (Days at sea) and the median number of bycaught dolphins per year. The two latter parameter may be the source of the variation between years of the number of bycaught dolphins. For example, the (posterior median) estimate is of > 600 dolphins in 2017 and it downs to < 100 in 2018. Total effort in DaS in the Bay of Biscay (divisions 8.a and 8.b) in the first ten weeks of 2017, when bycatch risk was highest, is twice the value of total effort in 2018. The median number of dolphins involved in a bycatch event is 2017 was also twice the number in 2018 (2 and 1 respectively, Table 5.1). The bycatch risk being equal for the ten first weeks of these two years, the estimate for 2017 is expected to be at least four times that of 2018.

Variations in the estimated number of bycaught dolphins through the model-based approach can be related to the estimated number of bycaught dolphins obtained from the strandings estimations. Stranding records are an independent source of data for estimating the number of bycaught dolphins (Peltier et al. 2016). Reverse drift modelling allows to infer the death location of each stranded dolphin showing bycatch evidence for each month between 1990 and 2020. Observed stranding numbers for each month can be corrected by both stranding and buoyancy probabilities (Peltier and Ridoux 2015). Reverse drift modelling cannot disaggregate estimates by métiers or fisheries but it provides an independent estimates of total mortality due to bycatch in the study area. The bycatch mortality estimated for the PTM should be lower than the total estimated from strandings. This was verified for all years save for 2010, but uncertainties were large and credibility intervals had a large overlap. By looking at the correlations between estimated number of bycaught dolphins for the stranding-based and the model-based approaches, it appeared that the correlations were significant between the two estimated series, especially after accounting for a lag due to drift (Figure 5.19). The magnitude of the correlation decreased between 2005 and 2020. One interpretation is that of a change in the relative contribution of PTM in total dolphin mortality over time, with PTM having a smaller impact on common dolphins in recent years compared to the 2000s.

## Within-Year variations

For each year, the bycatch risk was estimated to be higher in the winter part of the years (i.e. week 1 to 10). This pattern is largely congruent with the pattern seen in strandings of common dolphins in the Bay of Biscay (Dars et al. 2020). A limitation of stranding data is that it does not allow to discriminate between fishing gears (Peltier and Ridoux 2015). The Obsmer data in contrast included

The mortality trend obtained by the model based approach seems to be rather constant while the stranding based trend seems to increase in recent years. This suggests that the increase in mortality estimated from the strandings is not due to PTMs.

The correlation between the two trends is stronger when considering a 2-week delay in stranding-based estimates. This correlation decreases over the period considered. The contribution of PTMs to total bycatch mortality may be lower than before.
geolocalized bycatch events, with a spatial resolution at the level of ICES divisions kept for analysis. For the PTMs flying the French flag, we could identify divisions 8.a and 8.b from week 1 to 10 for nearly each year as the ones with the highest risk of bycatch. For each year, the number of bycaught dolphins can also be estimated at the ICES division level for each year by the ICES Working Group on Bycatch (WGBYC). Using data collected by onboard observers collected between 2005 and 2017, bycatch rates for ICES divisions on the continental shelf of the Bay of Biscay were estimated with ratio-estimators (ICES 2019). These estimates are not produced at the week-level, but (ICES 2019) also identified divisions 8.a and 8.b as the ones with the highest bycatch in mid-water trawls for common dolphins over the period 2005-2017 (Page 120). ICES (2019) estimated yearly rates ranging between 0.285 to 0.372 dolphins per DaS and warned against extrapolation given the low observer coverage). Our model-based approach overcomes this limitation (Authier, Rouby, and Macleod 2021) and allowed to identify, within each year, that weeks 2 to 6 were the ones with the highest bycatch numbers for both divisions 8.a and 8.b. This results were concomitant with the seasonal stranding pattern observed each year on the French seashore (that is, winter strandings; Peltier et al. 2013): around 80% of all common dolphin strandings on the French Atlantic seashore is observed between the end of January and the beginning of March.

A key feature of our model-based approach is how it leverages correlations between bycatch risk, haul duration and number of hauls per DaS (Figure 5.16). Some of the correlations are expected, such as the negative correlation between haul duration and the number of hauls per DaS. If more hauls are performed in a Day at Sea, so they are to be shorter and vice versa. However, average haul duration is not constant within a year, with the variations reflecting change in the commercial fish species targeted by PTM at different time of the year. These variations at the week-level were quite substantial, and were taken into account when estimating bycatch in our model. There was a positive correlation between haul duration and bycatch risk, with at least a two-fold increase in the later when haul duration exceeds 5 hours (Figure 5.16). This was particularly evident in weeks 1 to 5 in 2017 and 2019 and week 48 in 2016. The latter was due to a single fishing trip with 5 hauls that lasted > 10 hours, each of which resulted in a bycatch event. We recommend, in light of the within-year pattern in haul duration (Figure 5.16), to investigate management actions and mitigation measures on limiting haul duration in winter to assess whether bycatch may also be reduced. The ability of our approach to propose more refined recommendations is a considerable advantage to reduce the impact of bycatch on the common dolphin population. However, it must be balanced against the fact that it has some

The current ICES mortality estimates do not make a fine spatial-temporal distinction. Our approach allows this distinction to be made for PTMs.

There seems to be a correlation between the fishing duration and the bycatch risk with a threshold effect. Beyond 5 hours of fishing, the bycatch risk is greatly increased. Restricting fishing duration may be a mitigating measure.
limitations.

## Benefits and limitations of the method developed and used

The case study presented in this chapter is based on the application of the Authier, Rouby, and Macleod (2021) approach. This approach was developed to deal with low quality on-board observer data. A crucial assumption associated with this approach is the availability of accurate fishing effort data to perform the post-stratification step. The measure of effort chosen was the DaS, according to the current practice of international working groups. (e.g. ICES 2019). This guided the joint modelling of haul-level risk, haul duration and average number of hauls per DaS. This modelling choice has proved successful for PTM but cannot be developed to the same extent for other métiers, in particular for passive gears such as gillnets and trammelnets. In the latter case, a better measure of towing effort is soak time which takes into account the length and height of the net. This information can be difficult to collect, making it difficult to obtain an accurate measure of effort. This makes post-stratification difficult, but it is essential for the application of the model-based approach. The model developed for PTM may not necessarily transfer seamlessly to other métiers.

Since 2021, PTM flying the French flag are required to use deterrent acoustic devices (pingers¹³). If these devices are efficient to reduce bycatch risk, this may be taken into account in the model, by adding a covariate in Equation 5.3. This requires post-stratification with the covariate using a fishing effort that explicitly takes the covariate into account. Currently, the model can implicitly account for the effect of pingers because it accounts for variations in bycatch risk from year to year and within a year. The large-scale deployment of pingers in 2021, if effective in reducing risk, will result in a lower estimated risk than in previous years. In other words, the model does not necessarily have to take into account all covariates at the haul level, as long as the objective is prediction rather than explanation (Authier, Rouby, and Macleod 2021). Explicit consideration of the pinger effect is only necessary to make sense of variations in risk between and within years, but not necessarily to estimate these variations.

While Authier, Rouby, and Macleod (2021) concluded on increased accuracy of using regularized multi-level regression with post-stratification to estimate bycatch with observer data, they also found that estimate precision was low. This was also the case in this study (Table 5.4). A simple way to increase precision is to include Electronic Remote Monitoring observations in Equations 5.7 and 5.8. Doing so would result in increased precision as a

<figure><img src="image_106.png" /><figcaption>Figure 5.21: Photo of Pinger used in Bay of Biscay.</figcaption></figure>

©Parc marin d'Iroise.

13: https://www.mer.gouv.fr/protecting-cetaceans-annick-giradin-present
greater number of hauls (and possibly DaS) would be monitored. Even if the estimate precision is low, it may provide a sufficient precision to inform management procedures algorithms conjointly to other data sources and biological indicators.

## Implication for conservation

The large additional mortality due to anthropogenic activities on North-East Atlantic common dolphin population triggered a dedicated working group on emergency measures in 2020: the workshop on fisheries emergencies measures to minimize bycatch of short-beaked common dolphins in the Bay of Biscay and harbour porpoises in the Baltic Sea (WKEMBYC) took place remotely in spring 2020 (ICES 2020b) and informed an ICES advice that same year (ICES 2020c). This advice led to an infringement procedure issued in July 2020 against France for failing its obligations under the Habitats Directive which lists the common dolphin as a species requiring full protection on its Annex IV. The same day, the Paris Administrative Court of Justice condemned the French government for failing to transpose and apply in a timely manner the disposi-tions of the Habitats Directive and of Technical Measures regulating fisheries¹⁴ (in French). Following the unprecedented number of strandings in 2017, a national working group with fishermen, their representatives, government officials, Non-Gouvernmental Orga-nizations and academics was initiated to address the bycatch issue (Peltier et al. 2021). One recommended action was to improve esti-mates of bycatch due to high-risk métiers, and to develop adequate methodologies to analyse data from non-representative samples (Authier, Rouby, and Macleod 2021). The present work reports on a case study on PTM and operating for a large part in the Bay of Biscay, and to a lesser extent in the Celtic seas. Model-based estimates can (i) provide information on pressures on common dolphins, as required by the MSFD (EU 2008/56); and (ii) take into account the ICES recommendation to develop estimation methods to make best use of data already collected to inform timely man-agement (ICES 2020c). In the next chapter, we will see how it is possible to combine both estimated number of bycaught dolphins and demographic indicators to inform conservation algorithm and provide adapted tolerable bycatch thresholds.

14: http://paris.tribunal-administratif.fr/content/download/172866/1715763/version/1/file/1901535.pdf
# Conservation Biology and Management Procedures

6

The North-East Atlantic population of common dolphins is affected by bycatch which threatens its viability. The status of the population is of concern in relation to the conservation objectives set by the Habitats Directive and the MSFD. The conservation objectives are the basis for establishing management strategies. These management strategies determine the thresholds beyond which the conservation objectives will not be met. In the case of the common dolphin population, these thresholds are those of bycatch, or removals. The challenge is therefore to define management objectives that correspond to bycatch limits in order to build management strategies that will allow the population to achieve the conservation objectives. ASCOBANS has proposed as an interim conservation objective "to restore and/or maintain stocks populations at 80% or more of carrying capacity". Based on this objective, OMMEG proposed the following quantitative conservation objective for marine mammals in European Waters: "a population should be able to recover to or be maintained at 80% of carrying capacity, with probability 0.8, within a 100-year period". Appropriate management strategies must therefore be defined to achieve this conservation objective for the common dolphin population. In this chapter we will determine the appropriate management objectives, i.e. the removal threshold. This study is part of a Management Strategy Evaluation approach that takes into account as much uncertainty as possible to provide relevant threshold. This approach involves the implementation of an operating model based on a Pella-Tomlinson population dynamics and a Leslie Matrix fed by estimated the vital rates in Chapter 4. A control rule based on the removal limit algorithm is applied to define the bycatch threshold at each time step during the 100 years simulated in the management strategy evaluation. A performance statistic is obtained at the end of the 100 simulated years. This statistical indicator, coupled with the use of real data, will ultimately allow to decide on the bycatch threshold to reach the OMMEG conservation objective.

Publication in prep.

This chapter is associated with one publication:
Rouby, E. and Genu, M. Setting Bycatch thresholds for common
dolphin through Management Strategy Evaluation. In prep.

6.1 Introduction: Conservation biology and Management procedure 130

6.2 Methodology: Management Strategy Evaluation (MSE) 134

6.3 Results 142

6.4 Discussion 146
<figure><img src="image_107.png" /><figcaption>**Figure 6.1:** Flowchart of the whole MSE and Application case approaches. MSE runs for 100 years and allows the selection of a quantile. The application case runs for 1 year using real historical data and allows to obtain a distribution of removals. The quantile selected in the previous MSE framework is applied on the distribution to select a removal value.</figcaption></figure>
## 6.1 Introduction: Conservation biology and Management procedure

In Chapter 5, we quantified the bycatch pressure exerted by PTM on the common dolphin population. In Chapter 4, we have seen that the bycatch pressure (any fishing *métier*) seems to decrease the survival within the population and to lead to a decline in its growth rate. Bycatch is actually one of the major threats to marine mammals (Avila, Kaschner, and Dormann 2018) particularly for small-sized cetaceans (Reeves, McClellan, and Werner 2013; Gray and Kennelly 2018; Brownell Jr et al. 2019; Rogan, A. J. Read, and Berggren 2021) and is currently impacting the common dolphin population (Peltier et al. 2016) to a point of non-viability (Mannocci et al. 2012). The conservation of marine mammal populations requires understanding and assessing the consequences of bycatch (Wade et al. 2021). Two examples of unsuccessful conservation cases associated to bycatch were mentioned in General introduction Page 7 with both the Baiji and the Vaquita and it is necessary to avoid similar situations for other populations. To do so, it is necessary to define conservation actions to promote the viability of marine mammal populations, in our case the common dolphin. Abundance indices are used to detect declines but once a decline is detected it is often too late to prevent the population from going extinct or it takes longer to recover to initial conservation status (Wade 1998; Gerrodette 1987; Cooke 1994; R. Williams, A. Hall, and Winship 2008; Authier et al. 2020). This is a reactive approach that is not precautionary (see Margin-Figure 6.2 for diagram; A. J. Read 2010). The development and use of pro-active conservation tools is of interest to prevent the extinction of the North-East Atlantic common dolphin population. An example of a pro-active conservation tool is the Potential Biological Removal¹ (PBR) developed as part of the implementation of the Marine Mammal Protection Act (MMPA) (Wade 1998).

The MMPA is one of the pioneering texts in terms of managing marine mammal populations from a conservation perspective² and it defines clear conservation objective (hereafter CO). The three CO retained by the MMPA are: 1) Maintain the population above its Maximum Net Productivity Level (MNPL); 2) Allow the population to recover close to its carrying capacity (K); 3) Allow the population to recover considering a growth rate close to its maximum possible (Wade 1998). MMPA defines clear quantitative CO and lays out management objectives and remedial measures to meet the CO. The quantitative CO for the MMPA is the following: "a population will remain at, or recover to, its MNPL (typically 50% of the population's carrying capacity), with 95% probability, within a 100-year period." In the case of European waters, there is

<figure><img src="image_108.png" /><figcaption>**Figure 6.2: Classical conservation approach.** The classical approach involves 5 steps (A. J. Read 2010). This approach is not proactive and may lead to viability problems for dolphin populations (Wade 1998).</figcaption></figure>

1. PBR = 0.5 × Rmax × Nmin × Fr with Rmax as the maximum theoretical or estimated productivity rate of the population, Nmin the minimum population estimate in number of animals and Fr the recovery factor.

2. For more information on the MMPA: https://www.fws.gov/international/laws-treaties-agreements/us-conservation-laws/marine-mammal-protection-act.html
a lack of both: 1) a legally-binding CO for marine mammals and 2) management objectives with respect to human-induced mortality (ICES 2020c; Rogan, A. J. Read, and Berggren 2021). No clear and quantitative CO are defined regarding the impact of removals³ on marine mammals population in European Waters. To address this, ASCOBANS defined a fixed percentage of the abundance of marine mammals population as a threshold to set the CO regarding bycatch in European waters (as mentioned in Chapter 2 Page 35). ASCOBANS passed two resolutions to set removal limit of marine mammal populations:

► 1. "Unacceptable interactions" are, in the short-term, a total anthropogenic removal above 1.7% of the best available abundance estimate (Resolution 3.3 passed in 2000).

► 2. Bycatch reduction have to be reduced to less than 1% of the best available population estimate as an intermediate precautionary objective (Resolution 3.3 and Resolution 5.5 passed in 2006).

This fixed percentage abundance approach have been used by the European Commission (EC) to assess the capacity of European member states to achieve the "Good Environmental Status" required by the MSFD. This fixed percentage approach has the advantage of being very simple in its use and easily understood. To apply it, at least one abundance survey covering the area associated with the population to be conserved and estimates of bycatch removals is required. It remains limited for two main reasons: 1) Inability to take into account life-history traits and 2) Does not take into account potential errors and biases (Winship 2009). Abundance estimates may not be synchronised with bycatch estimates, limiting proactive conservation initiatives. A more conservative approach is to define explicit bycatch thresholds as management objectives in relation to quantitative CO. The OSPAR Marine Mammal Expert Group (OMMEG) provided quantitative interpretation of the ASCOBANS interim objective "to restore and/or maintain stocks/populations to 80% or more of the carrying capacity" (ASCOBANS 2000). The quantitative CO defined by OMMEG in European waters in the context of the MSFD is the following: "a population should be able to recover to or be maintained at 80% of carrying capacity, with probability 0.8, within a 100-year period.". Explicit bycatch thresholds can be produced in relation to this quantitative CO.

A management strategy is an agreed-upon set of rules for determining thresholds beyond which a CO runs the risk of not being met with unacceptably high probability (Winship 2009; Punt 2006; Bunnefeld, Hoshino, and Milner-Gulland 2011; I. C. Kaplan et al. 2021). Management strategy defines management objectives in the form of thresholds. Evaluating the relevance of management strategies (Management Strategy Evaluation, MSE) to

3: Removals are bycaught dolphins removed from the population.

In contrast to the USA, for which the MMPA defines a clear CO, Europe does not have a defined quantitative CO despite the ASCOBANS resolutions. OMMEG has provided a quantitative interpretation of the objectives set by ASCOBANS.

Conservation and management objectives should not be confused. A conservation objective is a goal to be achieved. It does not change. A management objective is a means to achieve the conservation objective. It is adaptable and can be reevaluated.
achieve COs is an issue that can be addressed through modelling and simulations (Cooke 1994; Hilborn and Mangel 2013). MSE requires generative models that generate (synthetic) data that are similar to observed data and it also requires real currently available data. Models associated with MSE are data-generating mechanisms. They can simulate a population dynamic taking into account removals effect (i.e. Bycatch effect) in addition to natural processes such as density dependence. These models allow scientists to assess the performance of management measures in hypothetical scenarios in order to set efficient management objectives. The management objectives will be assessed using real observable and available data (i.e. abundance and bycatch data with their uncertainties). Taking into account uncertainty at several levels (e.g. underlying model, potential biases and observed data) is of importance to measure the robustness of management strategies and choose the better one. A good management strategy should cover all aspects of management in accordance with pre-specified objectives. The thresholds regarding conservation of marine mammals are removal limits. A removal limit is an annual maximum number of animals whose removal would not result in excessive depletion of the population. MSE requires several components, including (list from Genu et al. (2021)):

1. One or Several Unambiguous quantitative CO.

2. A data simulator (or operating model) to emulate the dynamics of the marine mammal population and the effects of anthropogenic activities on this population.

3. A control rule, whose computation accounts for the expected quantity and quality of observable data, to set a removals limit beyond which the impact of human activities runs the risk of failing the aforementioned CO.

4. Performance metrics (or performance statistics), necessarily context-dependent and reflecting the trade-off between the potentially multiple CO defined previously.

These elements are necessary to project forward in time the pop-
ulation dynamics (i.e. the number of animals at each time step
according to population models operating within the data simu-
lator). For each management strategy the control rule is applied.
Finally the performance metrics are monitored and assessed with
respect to the CO. Items 1. and 4. should be agreed by all stakehold-
ers and scientists should not be expected to set the CO (Mangel et al.
1996). Items 2. and 3. are more relevant to scientists whose task is
to test a wide range of realistic scenarios in order to protect the
management strategy against uncertainties and potential biases in
the available data. MSE involves a tune *via* simulations. Tuning in
the MSE context means to: find, with a large number of simulations,
parameter values of the control rule that meet the CO (Genu et al.

The application of a management ob-
jective (e.g. removal limit) is a man-
agement strategy. The management
objective is implemented through a
strategy. To assess the adequacy of
strategies, the MSE framework can
be used. This framework involves the
definition of an operational model
and a control rule.
2021).

There is a need to produce sustainable bycatch thresholds for the common dolphin population through an approach that takes into account uncertainty by involving several types of data. This approach should be part of a management strategy evaluation process to enable the chosen bycatch limit to meet the conservation objectives set by OMMEG in European waters. In this chapter we will see the MSE procedure applied to the case of the common dolphin in the Bay of Biscay. This application is made in order to define bycatch thresholds for this population calculated according to the achievement of the quantitative CO set by OMMEG from the work of ASCOBANS. The work carried out in demography on the common dolphin (Chapter 4) makes it possible to inform the population dynamics model constituting the operating model of the procedure. The recent development by Genu et al. (2021) allows for the calculation and redefinition of management objectives during the process. In line with the steps outlined above, we will start by choosing the CO for this population (OMMEG CO). Then we will see the implementation of the operating model (Pella-Tomlinson, General introduction Page 12). Then we will detail the control rule (Removal Limit Algorithm, RLA) and finally the choice of a performance statistic (Quantile that define removals). These steps are described in the baseline scenario case (no bias in abundance and bycatch estimates). We will then explain the robustness trials that aim to measure the effectiveness of the management strategies in achieving the CO when bias is introduced (bias in abundance or bycatch estimates). After selecting the most appropriate performance statistic for the common dolphin through simulations, we will use this statistic to define a bycatch threshold. This definition is based on the use of real data in the procedure. This work was carried out in close collaboration with Mathieu Genu whose role was to code all the analyses as part of the development of the RLA package to implement the procedure described in this chapter. This work is a preliminary analysis and does not yet aim to give a limit. The work was carried out taking into account the Bay of Biscay which does not correspond to the extent of the management unit recognised for the common dolphin population. This choice was made to allow the use of current data from the Bay of Biscay and to start exploratory analyses. The following sections present the method and preliminary results.

We used an MSE framework to conduct exploratory analyses to define a removal limit (i.e. management objective) that meets the conservation objective set by OMMEG.
## 6.2 Methodology: Management Strategy Evaluation (MSE)

### Demographic adjustment to make the population stationary

The vital rates obtained in Chapter 4 do not allow the population to be viable (Figure 6.4 with constant effect = 0). To enable it to be viable an adjustment on survival was made through the adding of a constant. Several values were tested from 0.1 to 2.0 and the one that allows the population to be stationary was 0.65 (Figure 6.4). This constant was then added in the joint modelling⁴ and allows the estimation of age-specific survival rates (see Figure 6.3) that allow the population to be stationary. Using age-specific survival rates without adjustment, the population becomes extinct. The adjustment is made to obtain a stationary population. The adjusted age-specific survival rates are assumed to be representative of the natural survival of the population.

<figure><img src="image_109.png" /></figure>

<figure><img src="image_110.png" /><figcaption>Figure 6.3: Age-specific survival rates obtained taking into account the constant effect on survival.</figcaption></figure>

4: The constant effect (c) on estimation is applied on the survival model from Chapter 4 with a logarithmic transform to be consistent with the modelling of age-at-death $x_{it}$.

$$ \log x_{it} = (\mu + \log c) + \sigma \times \epsilon_i - \frac{Z_i}{\beta} \quad (6.1) $$

Figure 6.4: Remaining population over 100 years with constant effect on survival. The constant effect was applied to found the value for which the population remains stationary. With a survival intercept to which is added a constant effect equal to 0.65, the population remains stationary over time. This effect is then used in the estimation of the survival model to obtain the age-specific survival rates of a stationary population which will be used to inform the Pella-Tomlinson population dynamics model.

### MSE base case scenario

#### Item 1. Choice of conservation objective

The conservation objective used here for the common dolphin population in the Bay of Biscay is that of OMMEG: "attain or maintain the population at an abundance equal to 80% of its carrying capacity (i.e. K) over 100 years with a probability of 80%"
(OMMEG 2021), to be in line with ASCOBANS Resolutions 3.3 and 5.5.

## Item 2. Operating model and Pella-Tomlinson population dynamics

The operating model is the functional part of the Management Strategy Evaluation framework. The operating model corresponds to the true-state of Nature, the state of our common dolphin population that changes every year. The population's state is defined through its abundance and its depletion level relative to its carrying capacity. Changes are done through the arrival of newborns, the birth rate, the sex ratio, the proportion of the population that is mature, the age-specific survival (and implicitly mortality) that contributes to the removal of individuals and the removals on the population due to bycatch. No migration is assumed and taken into account.

The population model is age-aggregated with a maximum lifespan for the common dolphin population set at 24 years. The population model is based on a Pella-Tomlinson density dependence process (Pella and Tomlinson 1969):

$$N_{t+1} = N_t + r_{\max} \times N_t \times \left(1 - \left(\frac{N_t}{K}\right)^{\theta}\right) \quad (6.2)$$

Here, the MNPL is determined by $\theta$. In the case of $\theta = 1$, population dynamics is given by a logistic growth curve with MNPL = $K/2$. The value of the birth rate at time $t$ ($b_t$) is density dependent. Birth rate is the mean number of offspring per female assuming a sex ratio of 1:1. For example, a birth rate of 0.5 signify that each female produce 0.5 newborns per year. So 1 newborn is produced every two years by each female. Birth rate at time $t$ is expressed as follows:

Parameters of Equation 6.2: $N_t$ as the abundance at time $t$, $r_{\max}$ as the maximum population growth rate, $K$ as the carrying capacity and $\theta$ as the shape parameter.

$$ \begin{aligned} \text{birth rate}_t &= b_k + (b_{\max} - b_k) \times \left(1 - \left(\frac{N_t}{K}\right)^{\theta}\right) \\ &= b_k + (b_{\max} - b_k) \times \left(1 - (D_t)^{\theta}\right) \end{aligned} \quad (6.3) $$

The number of newborns in the population is the sum of the product of the population birth rate $b_t$, the age-specific abundance ($N_{at}$) calculated with the age-specific survival and the age-specific proportion of mature females $M_{at}$ in the population.

$$ \text{newborns}_t = \sum b_t \times N_{at} \times M_{at} \quad (6.4) $$

Birth rate ($b_t$) calculation on Equation 6.3 depends on the birth rate at carrying capacity $b_k$, the maximum birth rate $b_{max}$, the real simulated population size at time $t$ ($N_t$) and the carrying capacity $K$. The depletion level of the population is assumed as: $N_t/K$ and corresponds to the abundance of the population expressed in percentage of $K$.
The population model updates the total abundance $N_t$ and the age-specific abundance $N_{a_t}$ within a time-step of 1 Year. Each year the processes on the population dynamic are applied as:

► 1) Anthropogenic mortality is applied through bycatch removals on each age-classes of the population.

► 2) Natural mortality is applied through the age-specific survival rate on the population.

► 3) Animals reproduce through the age-specific fecundity rate and the abundance of age-class 0 is populated.

► 4) Increment the time by one year.

Processes on the population are shown in Figure 6.5 in the operating model block. The conceptual diagram of population dynamics from the demography chapter (Chapter 4) is taken up and adapted for the current study.

Since the initial depletion level of the population is unknown, 6 initial depletion levels were simulated simultaneously. The aim is to cover the population dynamics given several initial depletion intervals⁵ The population model runs for 6 years. After 6 years, the first control rule (Removal Limit Algorithm) is applied.

5: Interval for initial depletion levels (% of K) are: (0.3, 0.4] ; (0.4, 0.5] ; (0.5, 0.6] ; (0.6, 0.7] ; (0.7, 0.8] ; (0.8, 0.9].

## Item 3. Control rule with the Removal Limit Algorithm

The control rule corresponds to the Removal Limit Algorithm (RLA) which is expressed as follows:

$$ 
\begin{aligned}
\text{RLA}_t &= \hat{N}_t \times r \times \max(0, D_t - \text{IPL}) \\
&= \hat{N}_t \times r \times \max(0, \left(\frac{N_t}{K}\right) - \text{IPL})
\end{aligned}
\qquad (6.5) $$

RLA_t is the outcome, the maximal tolerable threshold removals that must be applied to the population in order to reach the CO. It is a management objective that set the removals applied each year on the operating model. It is reassessed each cycle of 6 years since every 6 years a survey campaign is simulated to fill in the value of $\hat{N}_t$, the estimated abundance of the population (in the base case scenario).

The estimation of RLA_t is carried out in a Bayesian framework. Its posterior distribution is obtained from the estimate of the posterior distribution of both $D_t$ and $r$ (Genu et al. 2021). From the posterior distribution of RLA_t a value must be chosen. It is the tuning part of the MSE. The tuning is a decision analysis. A quantile value within the posterior distribution is chosen in accordance with the achievement of the CO. A diagram illustrating the effect of quantile choice is shown in Figure 6.6. The choice of quantile allows

Parameters of the Equation 6.5. RLA_t is the bycatch threshold at time t, r is the population growth rate, $D_t$ is the depletion level of the population at time t and the IPL is the Internal Protection Level.
Management Strategy Evaluation

<figure><img src="image_111.png" /><figcaption>**Figure 6.5:** Diagram of the Management Strategy Evaluation associated with common dolphin population. The diagram shows both process (here blocks) that make up the Management Strategy Evaluation. The first is the operating model. This process involves a population dynamic that runs every year for 100 years. The second is monitoring in relation to the definition of the management strategy. Abundance monitoring is carried out every 6 years and bycatch monitoring every year. For each 6-year cycle a new bycatch limit is defined.</figcaption></figure>
to set management objectives: the tolerable bycatch threshold is
calculated by choosing a quantile from the Equation 6.5.

<figure><img src="image_112.png" /></figure>

**Figure 6.6: Diagram of tuning.** The frequency of removals thresholds is depicted as a posterior distribution. The choice of quantile affects the definition of the management objective. Given the posterior distribution of RLA<sub>t</sub>, a quantile value of 30% involves a tolerable bycatch threshold of 500 individuals per year. A quantile value of 70% involves a tolerable bycatch threshold of 2500 individuals per year. The choice of the quantile value therefore defines the management objective. This choice is made in order to achieve the CO.

One important element is the presence of the Internal Protec-
tion Limit (IPL). The IPL is set to 0.54 = 54% of K. If the abundance
of the population at time t, N<sub>t</sub> < IPL, so RLA<sub>t</sub> = 0. This protection
is set in order to limit removals to 0 if the population is too much
depleted. One essential assumption of this part is that the removal
limit is always respected (in the base case scenario) for the 6 next
years.

The steps constituting the items of the operating model can
be completed with those of the control rule:

**Operating model:**

► 1) Anthropogenic mortality is applied through bycatch re-
movals on each age-class of the population.

► 2) Natural mortality is applied through the age-specific survival rate on the population.

► 3) Animals reproduce through the age-specific fecundity rate and the abundance of age-class 0 is populated.

► 4) Increment the time by one year.

**Control rule:**

► 5) From the estimated abundance and simulations estimate the distribution of possible removals.

► 6) From this distribution, determine the value corresponding to the quantile applied.

**Item 4. Performance Statistic, Quantile for setting Removals**

Originally, the Management Objective set by ASCOBANS
was translated by the OSPAR Expert Group on Marine Mammals
(OMMEG) into the following: After 100years, the depletion level
of the population should be equal to 80% of its carrying capacity
with a probability of 80%. At the end of the simulated 100 years, the depletion level of the population is established. If it is above 80%, the management strategy is assessed as viable. But it is necessary to choose one. The choice is performed on quantiles. The management strategy which is retained to be applied is the one associated to the maximum value of quantile that allow to reach the management objective. For example, if quantile 0.2 is associated to a depletion of 90% after 100 years and the quantile 0.3 to a depletion of 82% after 100 years, the retained quantile will be 0.3. Because it allows to reach the management objective set by OMMEG under ASCOBANS. The selected quantile can then be applied in the Management Strategy Application as it will be depicted latter. Figure 6.1 summarises all the steps needed to apply the MSE for the base case scenario which is the case where we suppose that we can fully trust our abundance and bycatch estimation data. However, it is not always necessarily the case. This is why we can take into account some bias in these estimations.

## Robustness trials

We considered an abundance estimate bias and a bycatch estimate bias in order to represent some cases when we can not trust the data (abundance and bycatch) that we have. The interest of using these bias scenarios in addition to the base case scenario is that it allows to be more cautious while choosing a quantile to perform the real case Application and finally obtain removal limit threshold.

The first bias scenarios are based on the abundance estimate. As depicted in Figure 6.1 and Diagram 6.5, an abundance estimate is done every six years (simulating a SCANS campaign) and from this abundance estimate, the new removal limit is set and revised six years later. We studied the case when abundance is underestimated by half (e.g. 120,000 dolphins estimated instead of 240,000) and the case when abundance is overestimated by half (e.g. 360,000 dolphins estimated instead of 240,000). We choose this values since they cover the confidence interval of SCANS II + CODA and SCANS III surveys⁶ which is approximately equal to 0.4.

The second bias scenario is based on an underestimate of the bycatch. The population model updates the abundance of each age-class depending partly on bycatch removals which are carried out on each age-class. If the management strategy is conducted with the total abundance obtained from the population dynamic model given a set of removals without bias (base case scenario), the removal limit set during the management strategy will be set according to the bycatch removals assumed as true. But it is possible that we underestimate these bycatch removals and

Two robustness trials were performed:

1) Over-estimation of dolphin estimates by a half

2) Underestimation if bycatch estimates by a factor of 3
These robustness trials are preliminary. Further trials need to be carried out. The values have been set to explore the impact of bias on the final result.

6: See Chapter 2, Page 36 for information about SCANS and CODA and see Supplementary Chapter 7.4 on page 201 for all abundance estimates.
set inappropriate removal limit. To cover this scenario we set an underestimation ratio of bycaught individuals equal to 3. In that scenario, the estimated number of dolphins bycaught per year is three times lower than the actual number.

Each of the two robustness trials considered here is intended to illustrate the difference in the choice of quantiles at the end of the analysis. They need to be complemented by other robustness trials to cover several possible situations. The definition of robustness trials is the last part of the MSE framework. Figure 6.1 shows the MSE steps given the base case scenario but the same scheme is applicable for the bias scenarios. After running the simulations over 100 years, given several types of scenarios, we move on to the final stage of this study. The interest of having carried out an MSE approach is that a quantile value has been defined which is adapted to the achievement of the Management objective and which respects a precautionary principle without being totally conservative.

## Application on common dolphins in the Bay of Biscay

The final step of the whole procedure is to apply the MSE approach using real historical abundance and bycatch data. The operating model is used once with the historical data. Then the control rule is applied to obtain the posterior distribution of possible removals. The performance statistic (quantile) selected from the tuning in MSE base case scenario and robustness trials is applied on the removals distributions to choose a bycatch threshold (see Real case Application on Figure 6.1).

The abundance data that are used are related to the Bay of Biscay. This scale is not representative of the common dolphin management unit. We made this choice in order to start exploring the analyses and setting up the approach. We also made this choice on the basis of the available data for two reasons. First the demographic rates that were obtained in the Chapter 4 were more linked to individuals that are associated to the Bay of Biscay than to the whole North-East Atlantic. It seems to me more accurate to use a similar spatial scale. Secondly, bycatch estimates are obtained from the analysis of drifting trajectories of stranded individuals. Here again, the individuals are mainly associated with the Bay of Biscay area. Thirdly, the management and conservation issues for the common dolphin in the North-East Atlantic are mainly linked to bycatch pressure in the Bay of Biscay. It is for these reasons that I believe it is currently more appropriate to use the estimated abundances in relation to the Bay of Biscay only. The tolerable bycatch thresholds that will be obtained subsequently will therefore be related to the Bay of Biscay "population". The abundance estimates

The practical application is based on the use of real historical data in addition to the demographic data used previously. The data used are abundance estimates and bycatch estimates.
used for the applied case are obtained from three surveys that represent two time points of abundance assessment separated by about 10 years: SCANS II + CODA and SCANS III⁷. Abundances were estimated from observation by plane and boat. Where Scans III directly combines the two, SCANS II does not cover the whole area and the CODA ship survey was needed to complete the estimates. Abundance estimates were made depending on areas composed by "blocks" Figure 6.7.

<figure><img src="image_113.png" /></figure>

The blocks were not the same for the SCANS II + CODA and SCANS III campaigns. Abundance⁸ data are provided by blocks. To allow for better cohesion we have recalculated the abundances taking into account the ICES divisions. We recalculated the abundances by multiplying the surface area of each ICES division by the estimated common dolphin density for the campaigns SCANS II + CODA and SCANS III (Figure 6.8).

<figure><img src="image_114.png" /></figure>

The other historical data that was used to perform the Real case Application was the historical bycatch estimate made from the strandings between 1990 and 2020. Theoretically, it is possible to use French pair trawlers (PTM) estimates that were produced in

7: ObSERVE is an aerial campaign conducted in the Ireland's EEZ whose abundance surveys are not useful for application in the Bay of Biscay. More information on ObSERVE on https://www.gov.ie/en/publication/12374-observe-programme/

**Figure 6.7: Survey platform and survey blocks of abundance surveys campaigns.** SCANS I is not shown here as it was not used to obtain the abundances. Between campaigns, the blocks sampled are different. A recal-culation of the abundances was nec-essary. This calculation was carried out within the ICES divisions. It took into account the density obtained in each block and the surface of the ICES division.

8: The abundance data used are not representative of the entire common dolphin management unit. The results obtained at the end will not be applicable to the management unit. This choice was made to facilitate the exploratory approach.

**Figure 6.8: Common dolphin abundance used for this study depending on the ICES divisions.** Abundances were recalculated using both ICES divisions specific surface and the estimated density in the SCANS II + CODA and SCANS III block.
Chapter 5. But the estimation would be made according to one of all the fishing métiers operating in the Bay of Biscay. Other fishing métiers than PTM could have an impact on common dolphin in the Bay of Biscay. Within the French fishing fleet, gillnetters (GTR) and gillnetters (GTR) may be associated with common dolphin mortality in ICES divisions 8.a and 8.b. The approach used in Chapter 5 is currently only applicable to PTM. In the meantime, it is possible to use the estimates made from the strandings. To estimate bycatch, only common dolphins with evidences of contact with fishing nets were considered (Kuiken 1994) as well as the common dolphins found stranded during multiple or unusual mortality events related to fisheries. These multiple events are defined when a high number of strandings occur in a restricted area with a common cause of death. The threshold is defined at 30 cetaceans over 10 consecutive days recorded along a maximal distance of 200km in the Bay of Biscay and 10 individuals per 10 days per 200km of coastline along the coast of the western English Channel (Peltier et al. 2014). From the stranding location of the dolphin, a deterministic drift trajectory is inferred by a reverse drift model (developed by Méteo France) to determine its likely area of death at sea. The area of death is included in a grid cell. The number of dead dolphins in each cell is then corrected by the cell-specific probability that a carcass from that cell washes ashore (Peltier et al. 2013). The number of animals estimated dead at sea is available in supplementary material.

Theoretically it is possible to use the estimates made in Chapter 5. For the time being they only concern one fishing métier and are therefore unusable.

## Baseline Management Strategy Evaluation assumptions

* The vital rates estimated in Chapter 4 are representative of those of the managed population.

* There is one closed and stationary population in the Bay of Biscay.

* The Conservation Objective produced by OMMEG is accepted by all parties and remains the same over the 100-year period considered.

* Abundance surveys are conducted every 6 years within 100 years.

* The Control rule used to produce the removals estimates remains the same over the 100-year period considered.

## 6.3 Results

Results are first presented for Management Strategy Evaluation comparing the depletion trajectories obtained under each quantile and each scenario. Then the population dynamic under
the real case application and the removals threshold obtained are presented.

# Management Strategy Evaluation

## Demographic dynamics of the base case scenario

Figure 6.9 presents the simulated estimated abundance and associated applied bycatch (limit set by RLA + stochastic variation) along time-line for the base case scenario depending on each quantile. The value of depletion is estimated every year for 100 years. Results were averaged for all initial depletion levels. Only quantiles 0.6 to 0.8 are displayed because under 0.6 the CO is always reached. For each quantile the estimated abundance increases before reaching a plateau. The plateau is not reached after the same number of years depending on the quantile. All abundance estimation profiles are above the management objective (80% of K) except for the quantile at 0.8. The quantile whose profile is closest to the management objective is the one at 0.7.

<figure><img src="image_115.png" /></figure>

**Figure 6.9:** Base case scenario abundance and removals over years associated with quantiles. Only quantiles 0.6, 0.7, and 0.8 are depicted here. These are the most interesting to look at because it is at this point that the 80% threshold defined for K is reached. All quantiles below 0.6 meet the objective after 100 years. Here it is the quantile at 0.7 that should be retained. Because it is the one that best fits the conservation objective. It is associated with annual removals that fluctuate around 4100 at the end of the 100-year simulation.

Removals are also displayed. Their value increases as the quantile increases with for example a maximum value of 3500 at the 0.6 quantile and a maximum value of 4500 at the 0.8 quantile. Likewise for the simulated estimated abundance, the profile shows an increase and a plateau which is reached more quickly when
the quantile is high. The estimated bycatch value is about 4000 dolphins with the quantile at 0.7 after 100 years.

## Selection of the best quantile and associated removals

<figure><img src="image_116.png" /></figure>

**Figure 6.10: Robustness trials depletion levels for each quantile.** When we start the simulation, we do not know the depletion state of the population with respect to K. To cover this uncertainty, 6 initial depletion levels are simulated for each scenario and quantile. In the base case, the quantile to be retained is 0.7. In the case of the underestimation of abundance by a factor of 2, the quantile to be retained is 0.5. In the case of an underestimation of the bycatch by a factor of 3, the quantile to be retained is 0.3. The absence of curves for the last quantiles of the robustness trials is due to numerical problems.

As the initial depletion state of the population is not known, several initial depletion states were simulated. For each scenario (base case, underestimated abundance, underestimated bycatch), six initial depletion states (from [0.3, 0.4] to [0.8, 0.9]) were simulated under seven quantiles (from 0.2 to 0.8). At the end of the 100-year time series, the performance of each management strategy is assessed by comparing the value of each depletion curve to the 80% carrying capacity limit (Figure 6.10). The best 80% depletion values for each initial state are presented over a period of 1 to 100 years. For the baseline scenario, the best quantile to select based on the simulation is 0.7. This quantile value allows each depletion curve to meet the management objective of 80% K. For the underestimated abundance scenario, the best quantile to select
is the quantile at 0.5. For the underestimated bycatch scenario, the best quantile to select is the 0.3 quantile.

For each combination of quantile and scenario, it is possible to calculate the value of removals as a percentage of K (Figure 6.11). The resulting distribution of removal for all time steps (independent of the initial simulated depletion state) is summarised. The higher the quantile, the higher the removal values. These removal values are purely indicative and do not constitute the final result of the approach.

<figure><img src="image_117.png" /></figure>

**Figure 6.11: Removals to be applied in percentage of K.** Each level of removals relatively to K are depicted for each scenario and quantile. Lower the quantile, lower the removals. The scenario with the highest removals is the base case. The scenario with the lowest removals is the one that includes a bycatch estimation bias.

## Application on common dolphin population

<figure><img src="image_118.png" /></figure>

**Figure 6.12: Final removals posterior distribution.** The posterior distribution obtained in the Real Application case is depicted here. The maximum removal value is about 12500. The minimum is about 0. The quantiles that were obtained previously in the MSE approach involving the robustness trials allows to select removal value in the distribution.

The previous exercise (management strategy evaluation) aimed to obtain the quantile values for each scenario to achieve the CO. The application case aims to produce a distribution of removals using the actual data. The quantiles obtained previously are used to define the removals to be applied from this distribution (Figure 6.12). They range from 0 to approximately 13000 given
the last abundance estimate. From this distribution (Figure 6.1 management application part to a representation of this step in the framework.) it is possible to choose values given the quantiles selected under the MSE. The quantiles obtained from the previous MSE approach are presented as dashed color lines (Figure 6.12) and associated removals are given in Table 6.1. The difference between quantile at 50% and 70% is about a factor of 2 (approximately 2000 and 4000). The quantile at 30% is associated to the lowest number of removals which is about 656 considering the last SCANS III abundance survey estimation of common dolphins in the Bay of Biscay.

**Table 6.1: Removals to be applied.** Here the removals that allows the population to reach the OMMEG CO are depicted. The higher removal value is the one for the base case with 4163 removals per year. The mid removal value is the one for the underestimated abundance with 2088 removals per year. The lower removal value is the one for the underestimated bycatch with 656 removals per year. In order to be conservative and to be cautious with the population viability, the lowest removal value should be kept. Here the yearly removals applied to the common dolphin population in the Bay of Biscay must not exceed 656 individuals per year in order to reach the CO set by OMMEG.

<table><thead><tr><td>Quantile</td><td>Removals</td></tr></thead><tbody><tr><td>30%</td><td>656</td></tr><tr><td>50%</td><td>2088</td></tr><tr><td>70%</td><td>4163</td></tr></tbody></table>

## 6.4 Discussion

The aim of the study is to define a tolerable removals (i.e. bycatch) threshold (management objective) to achieve the conservation objective defined by OMMEG based on the ASCOBANS recommendations: allow the population to reach or maintain an abundance level equal to 80% of its carrying capacity (i.e. K) over 100 years with a probability of 80%. The implementation of this management objective will then constitute the population management strategy until the next abundance surveys.

For this purpose, demographic data were used considering an adjustment on survival to obtain a stationary population. Then the MSE procedure was applied involving the selection of the CO, and the combined use of the operating model and the control rule over a period of 100 years to obtain the performance statistic. Actual abundance and bycatch estimates were then used to compute tolerable removal thresholds for the population given the latest abundance estimates. The choice between these thresholds then depends on the quantile selected from the MSE approach and the robustness trials associated with the base case scenario.
# Conservation objective, operating model and control rule

## The choice of an appropriate CO

The definition of a clear and justified CO conditions the subsequent MSE approach. The CO is relative to K, the carrying capacity of the environment for the population. In practice it is impossible to know the exact value of this carrying capacity. But knowing its exact value is not limiting. The order of magnitude of the carrying capacity, however, is much more important. In the case of the common dolphin population, it is more important to decide between a carrying capacity in the tens, hundreds of thousands or millions. It is at this level that the maximum abundance surveys inform the CO. Defining the percentage threshold of this carrying capacity to be achieved is a more general management consideration. This definition directly determines the number of annual tolerable bycatches in the population. The higher the percentage, the lower the number of removals and *vice versa*. If the percentage set is very high, *e.g.* 95% of K, then only the lowest quantile values will achieve this CO (Figure 6.10). In this case, the tolerable number of bycatches will be close to 0. It may be difficult to implement this threshold in view of the socio-economic stakes of fishing in the Bay of Biscay. If the percentage set is low, *e.g.* 50% of K, then the higher quantile values will also achieve the CO. In this case the number of bycatches may be very high. The management objective will do little to limit fishing but may be dangerous to the viability of the population. It is possible that the population will be too depleted before the next abundance survey (6-year cycle) to redefine the management objective. The population could also be subject to an Allee effect that would cause its extinction beyond a given abundance threshold. It could be interesting to implement an Allee effect in the operating model to take this into account. This Allee effect could be expressed as a probability of encounter between males and females as an example. The question of how to define an appropriate CO that takes into account the conservation of populations and the management of human activities is a key issue at the border of several disciplines.

The exact value of K has little influence on the final estimate. It is the order of magnitude of K that matters. The choice of an appropriate CB must take this order of magnitude into account.

## Considerations on the Operating model

The MSE used here is based on multiple steps and fundamental assumptions. The operating model is based on the Leslie population matrix model with a Pella-Tomlinson density dependence for birth rates (Genu et al. 2021). The age-specific survival rates used to fuel this model are those obtain for the common dolphin population in the Chapter 4 adjusted to allow the population to be stationary. This construction makes it possible to obtain
a population that does not decrease without the application of bycatch. This adjustment is a strong assumption that does not go hand in hand with what is obtained in terms of demography on this population. The Leslie matrix population model used here does not admit explicitly immigration. This assumption is very restrictive since the common dolphin population in the Bay of Biscay seems not closed as the abundance surveys and demography suggest.

## Control rule and importance of data provision

The main interest of this work is that it is the first to be done with specific and locally determined demographic parameters on this population. This age-specific structured model allows the use of the RLA which is a possible control rule among others. Others control rules that can be used in the MSE are the PBR and the mPBR⁹. They did not require age-structured models but only the abundance of the population. The PBR can be calculated as it follows:

$$PBR = 0.5 \times R_{max} \times N_{min} \times F_{r} \quad (6.6)$$

From this control rule, a total of 4926 individuals per year were estimated as the limit of removals for the population in the whole North-East Atlantic area (ICES 2020b). The number of estimated bycaught dolphins for recent years is higher than this threshold. This threshold is close to the one calculated here, in the base case scenario, using the RLA for the Bay of Biscay. The calculation of the PBR relies on the management objective of the population reaching 50% of carrying capacity. This CO is the one implemented by the MMPA. It does not correspond to the CO proposed by the OMMEG group using the interim objective from ASCOBANS. The mPBR is suited for the ASCOBANS management objective and is associated to much lower bycatch thresholds with a removal limit around 900 individuals. This limit is closer to the estimation done using the RLA considering the scenario of a bycatch underestimation.

The use of RLA is less conservative when the accuracy of the data is higher and when there is more data. The removal limit determined when using the RLA on the Bay of Biscay population is close to that obtained from the PBR in the North East Atlantic. In this case, the use of population-specific demographic data allows for less conservatism. The use of the RLA in the framework of robustness trials also shows this particularly in the case of bycatch. In the base case scenario, bycatch estimates are assumed to be accurate. This may look like transparent bycatch declarations from fisheries. In this case, the tolerable removals on the population

9: mPBR is the modified PBR. It allows the estimation of tolerable removals thresholds taking into account the European conservation objective.

Parameters of PBR: $R_{max}$ as the maximum theoretical or estimated productivity rate of the population, $N_{min}$ the minimum population estimate in number of animals and $F_r$ the recovery factor.
of the Bay of Biscay is 4163 dolphins per year. In the case where the number of catches is underestimated by 3, this limit is greatly reduced by a factor of almost 7 (656 tolerable removals each year in Bay of Biscay). Robustness trials in that case shows that knowledge of bycatches makes it possible to be less conservative. It is therefore in the interest of fisheries to communicate transparently on the number of bycatches made per year.

## Application of the Management Strategy

In the case of the common dolphin in the Bay of Biscay, the most efficient management strategy to promote the sustainability of the population (without completely closing the fisheries) seems to be to set a limit of 656 bycatches per year. This threshold is associated with an underestimate of bycaught dolphins by a factor of 3. This scenario is extreme and other and less extreme robustness trials should be tested in the future. Within this study, this management objective seems to be the most precautionary given the scenarios tested and the current data. The base case scenario with a bycatch threshold of 4163 implies full confidence (i.e. unbiased) in the abundance estimates, bycatch estimates and in the demographic data in the operating model.

Bycatch estimates are currently made by modelling the reverse drift of carcasses (with bycatch evidence) found stranded on the seashore. It is difficult to be sure that these estimates are accurate. For example, they may be overestimated. The overestimation may be related to a possible misinterpretation of some bycatch evidence in favour of an over-representation of bycaughts individuals¹⁰. This hypothesis seems unlikely given that the evidence of incidental captures is well characterised for common dolphins (see Figure 2.19 in Chapter 2). In France, the examination protocol is well defined and applied by the members of the national stranding network. But the scenario of overestimation of catches will be tested. The scenario of underestimated catches foresees a maximum 656 catch limit. If a catch limit were to be applied, it would be 656 individuals per year. This is the lowest limit associated with the robustness trials. It is the most conservative value. Taking this threshold into account, the estimates we made in Chapter 5 are still below the threshold. Except if we take into account the 80% limit of the Bayesian uncertainty interval. But our estimates are for a single fishing métiers (Table 6.2).

Bycatch estimates are above this threshold in the Bay of Biscay, taking into account baseline stranding estimates. The current level of mortality reflects the current management strategy for the common dolphin population. This strategy involves the economics of the fishery, the viability of the common dolphin population, the

<figure><img src="image_119.png" /><figcaption>Figure 6.13: Cause of death for 2020.</figcaption></figure>

<footnote>10: As a reminder, 80% of individuals are considered dead due to capture in fishing nets on the basis of external examinations in 2020 and 82% on the basis of internal examinations.</footnote>
<table><caption>**Table 6.2:** Difference between the threshold of 656 removals and estimated bycatch of dolphins. Model based estimates correspond to the Chapter 5 estimates with the 80 Bayesian Credibility Interval. Cells are depicted in red if the difference is in favour of an excess of the bycatch threshold. Model based estimates are made taking into account PTMs.</caption><thead><tr><th>Year</th><th colspan="2">Model based (PTMs)</th></tr></thead><tbody><tr><td>2004</td><td>-656</td><td>-429 +478</td></tr><tr><td>2005</td><td>-656</td><td>-363 +761</td></tr><tr><td>2006</td><td>-656</td><td>-370 +647</td></tr><tr><td>2007</td><td>-640</td><td>-579 -437</td></tr><tr><td>2008</td><td>-644</td><td>-591 -466</td></tr><tr><td>2009</td><td>-473</td><td>-244 +164</td></tr><tr><td>2010</td><td>-652</td><td>-424 +338</td></tr><tr><td>2011</td><td>-647</td><td>-465 -21</td></tr><tr><td>2012</td><td>-633</td><td>-290 +534</td></tr><tr><td>2013</td><td>-643</td><td>-230 +890</td></tr><tr><td>2014</td><td>-656</td><td>-572 +272</td></tr><tr><td>2015</td><td>-654</td><td>-562 -210</td></tr><tr><td>2016</td><td>-601</td><td>-386 -273</td></tr><tr><td>2017</td><td>-500</td><td>-38 +759</td></tr><tr><td>2018</td><td>-655</td><td>-621 -493</td></tr><tr><td>2019</td><td>-597</td><td>-440 -215</td></tr><tr><td>2020</td><td>-652</td><td>-599 -466</td></tr></tbody></table>

<figure><img src="image_120.png" /><figcaption>**Figure 6.14:** Diagram of qualitative management strategy. Rather than consider a quantitative management strategy, strategies may be defined taking into account several socio-economic variables. A. Strategy which is more concerned with the sustainability of the fishing economy. B. Strategy which is more concerned with the viability of dolphin and fish populations.</figcaption></figure>

viability of the target fish species, public opinion and the certainty of estimates. In addition to these components, management effectiveness is also taken into account. A qualitative representation of the performance of these components for each management strategy can be imagined (Figure 6.14). This representation is adapted from (Punt et al. 2016). The current management strategy in the Bay of Biscay looks more like A. Despite ICES recommendations to implement time-area closures of fisheries, the government has chosen to favour the fishing economy. Public opinion is negatively influenced by the reports and the action of NGOs on the subject. If bycatch thresholds are defined and respected, then the management strategy could move from format A to format B. The sustainability of the dolphin population would be better taken into account, to the detriment of the fishing economy and in favour of public opinion.

This chapter closes the analytical part of the manuscript. At this stage we have explored different aspects:
* Chapter 2 on page 26 : Description of the common dolphin, its population in the North-East Atlantic, the associated conservation issues and the monitoring carried out on this population.

* Chapter 3 on page 43 : Development of a new survival analysis approach that takes into account covariate and random effects to estimate vital rates.

* Chapter 4 on page 70 : Estimation of vital rates for the common dolphin population, within the framework of the MSFD, which are estimated by taking into account various covariate and random effects.

* Chapter 5 on page 102 : Explanation and quantification of the bycatch phenomenon involving PTMs on a weekly, annual and ICES divisional basis.

* Chapter 6 on page 128 : Establishment of a tolerable by-catch threshold given the population's demographic viability, abundance, bycatch pressure and the CO set by OMMEG in European waters.

The following part is a discussion presenting three aspects related to this thesis. The first section deals with the general approach of the project. The second section deals with the North-East Atlantic common dolphin population. The third section deals with the cross-sectional monitoring.
DISCUSSION
<figure><img src="image_121.png" /></figure>
## 7.1 Project approaches

### Aspects related to the demographic study

This thesis project was carried out using demographic approaches. The first approach was to develop a new methodology for estimating survival, based on the Reed model, from age-at-death data. This approach is a useful contribution to the literature. Its use implies the inclusion of an inter-individual frailty parameter. In addition to allowing the modelling of a bathtub shape hazard curve, the use of the frailty parameter allows heterogeneity within the study population to be taken into account. It is difficult to interpret this parameter on its own and it is more accurate to consider it as a statistical device for robust estimation. However it may be possible to strengthen its biological interpretation. It would be interesting to consider two inter-individual frailty parameters. One associated with the calf component of the life history, the other with the non-calf or juvenile-adult component. This would allow a better representation of survival and mortality in the population, taking into account the different frailties at different periods of life. The survival estimation method is based on the use of age-at-death data and it is difficult to know whether the estimated vital rates (e.g. survivorship) are representative of those in the population. A possible solution would be to have both longitudinal monitoring data and age-at-death data for a population. Comparison of the longitudinal survivorship profile with that obtained through the cross-sectional survival approach would highlight the accuracy of the survival estimates. A limitation of this comparison is that, if there is a difference, the source of the difference cannot be defined. Is it a difference due to sampling or a difference due to the statistical approach or a combination of both? One way to compare would be to simulate longitudinal data as well as age-at-death data. Sampling and bias control would be provided by the simulation study. Any differences would be due solely to the statistical estimation of survival.

To control for sampling bias, stratified random sampling was useful. The use of this sampling framework allowed for greater confidence in the estimates of vital rates. When sampling in this way, the survival profiles estimated by the Siler and Reed models were almost identical. But this sampling scheme can only be used if

7.1 Project approaches ... 153

7.2 Definition of the North-East Atlantic common dolphin population ... 158

7.3 Population dynamics of elusive species ... 167

7.4 Conclusion and perspectives ... 175

A possible improvement in the estimation of survival is to take into account two frailty parameters. One parameter specific to the calf component and the other specific to the non-calf component. Such an improvement could be made through a mixture modelling approach.

To check the plausibility of the survival rates estimated by the statistical approach, it would be interesting to simulate both longitudinal and age-at-death data and compare the survival estimates.
there is a highly structured, centralised and coast-wide stranding Network. The French national stranding network is one of the most stable and unified stranding network in Europe. In Spain there is no single, centralised stranding network. Spain has many regional strandings networks along the mainland coast as well as separate networks for each of the Canary Islands. Germany has two stranding networks, while the Netherlands has four. Few European countries have a single stranding network. For those countries where this is the case, the network may not always be stable or long established (Figure 7.1. The French national stranding network has the advantage of being very reliable since the 1990s in recording all strandings on the Atlantic coast (Authier et al. 2014). Of the 25 stranding networks operating in North-Atlantic European waters, 11 routinely determine the reproductive status of females and 7 on *ad hoc* basis. With regard to age determination, data acquisition is variable since only 3 networks determine age as part of a routine analysis, while 13 do so on an *ad hoc* basis (see Margin-Table 7.1). It may therefore be difficult to generalise this sampling approach to the whole of the European level.

<table><caption>**Table 7.1: Frequency of analysis from stranded animals performed by European stranding networks.**
Table from ICES (2021a).</caption><thead><tr><th>Analysis</th><th>Routine</th><th>Ad hoc</th><th>Never</th><th>No reply</th><th>Total</th></tr></thead><tbody><tr><td>Female reproductive status</td><td>11</td><td>7</td><td>1</td><td>6</td><td>25</td></tr><tr><td>Maturity state</td><td>10</td><td>7</td><td>1</td><td>7</td><td>25</td></tr><tr><td>Diet</td><td>9</td><td>10</td><td>0</td><td>6</td><td>25</td></tr><tr><td>Male reproductive status</td><td>7</td><td>10</td><td>1</td><td>7</td><td>25</td></tr><tr><td>Age determination</td><td>3</td><td>13</td><td>3</td><td>6</td><td>25</td></tr></tbody></table>

<figure><img src="image_122.png" /><figcaption>**Figure 7.1: Period of activity of the networks that provided the starting date of their activity and the perceived stability of effort over time (as indicated by the respondents of the WGMME questionnaire).** Very few European countries exhibit a unique and long-term stranding network. Island, United-Kingdom, Belgium and France have the most structured stranding networks.

Figure from ICES (2021a).</figcaption></figure>

The possibility of including covariate and random effects in the analysis can provide a demographic indicator for the MSFD. At present, the criteria associated with obtaining demographic data for the MSFD descriptors are considered secondary. In other
words, Member States' efforts to fulfil the descriptors should not focus on these criteria as a priority. Criteria associated with obtaining abundance data, for example, are considered primary criteria that can provide information on the viability of marine mammal populations (including the common dolphin). In Chapter 4, we have seen that the viability of the population seems to have deteriorated over the period of time considered (1997 to 2019). The use of random effects taking into account the cohort effect allowed us to identify possible differences in survivorship and age at sexual maturity for females. These differences could be used to develop an operational demographic indicator to immediately inform the demographic criteria of the MSFD descriptors. This indicator could take the form of a percentage difference in vital rates between cohorts (Figure 7.2). Differences from a baseline or mean profile could be reassessed every 6 years in each MSFD cycle.

Estimating vital rates taking into account random effects can be a new indicator of demographic status within the MSFD. This would help to identify the possibility of depletion before it occurs.

<figure><img src="image_123.png" /><figcaption>**Figure 7.2:** Relative percentage difference of cohort vital rates around the mean profile from 1997-2019. Relative percentage difference in survivorship. From 2012, the survivorship of adults seems reduced compared to the mean survivorship over the period.</figcaption></figure>

The use of covariates also revealed a possible effect of bycatch on vital rates. This effect is quantifiable and could also be a pressure indicator for the demography of the common dolphin population to inform the MSFD. The effect of various anthropogenic pressures, such as pollutants, can be quantified at the demographic level for this population. For these reasons, I think that the criteria for assessing the "Good Environmental Status" of European Waters based on demography should not be considered secondary. Vital rates are important for knowing the viability of the population and it is now possible to quantify variations in the viability of the population. Vital rates provide earlier information on whether the viability of a population is impaired than abundance surveys. It is possible to identify a risk of depletion before it is observed. These criteria should be reviewed as primary criteria.

This demographic indicator could be used to fulfil the MSFD at the

Demographic criteria are considered secondary. They should be reviewed as primary criteria. Demographic information for elusive marine mammal species can provide more detailed information on impacts and viability.
sub-regional scale. The MSFD considers the evaluation of water quality on a regional and sub-regional scale. For the common dolphin, the regional scale is the "North-East Atlantic", while the sub-regional scale is the "Bay of Biscay and Iberian Coast" and "Celtic Seas". Using the individuals stranded on the shores of these two sub-regions, an operational indicator could be produced to inform the demography of the common dolphin within the two sub-regions.

## Aspects related to the conservation approach and management issues

This thesis project was carried out with the implementation of a management approach. The management approach is based in part on the use of bycatch estimates. The approach used in Chapter 5 to quantify bycatch from on-board observer data shows promise for using biased data. This approach can be used to add value to data that already exists but are considered biased. The ability to work with data that are considered biased should not prevent the collection of more representative data. Although promising for taking observation bias into account, the use of this approach is currently limited to PTMs. To extend its use, it would be necessary to be able to effectively calculate the fishing effort of passive *métiers* such as gillnetters. The use of this approach in the context of implementing management strategies currently seems to be limited since it cannot be generalised to all fishing *métiers*. On the other hand, it is interesting in terms of improving data that are considered to be of little use by the research community and providing information on the phenomenon of bycatch at a fine scale. The use of the joint model also makes it possible to highlight possible correlations or links between the various estimated variables (e.g. fishing duration and bycatch risk).

The removal threshold produce through the MSE approach does not include any socio-economic index as it is done for some fisheries evaluation models such as the Atlantis model. The Atlantis model includes all the components of marine ecosystems from physics to human exploitation (Audzijonyte et al. 2019; Figure 7.3). The use of these complex models may be limited by their complexity. It can be difficult to define the links between all parts of the model and to understand how these parts interact. Analyses in the context of an MSE incorporating RLA are already a complex task involving many assumptions. For this reason, I think that the development of more complex models may be interesting from a theoretical point of view but not necessarily interesting from an applied point of view. Discussions between stakeholders (fisheries, scientists and decision-makers) are already difficult with regard to the bycatch and removal threshold. In my opinion, we need to keep it simple

The estimation of the number of animals bycaught and associated variables makes it possible to valorise biased observer data. These data provide fine-scale spatio-temporal estimates. However, the use of this approach is limited by the availability of an accurate fishing effort dataset and sufficient bycatch observation data.

<figure><img src="image_124.png" /><figcaption>**Figure 7.3: Diagram of Atlantis model.** Models blocks are depicted such as relation between blocks. The model takes into account the whole ecosystem from the oceanography to the top of marine foodweb and human activities.</figcaption></figure>

Figure from https://niwa.co.nz/ecosystem-modelling-at-niwa/atlantis-ecosystem-model.
in order to facilitate understanding and decision making on the strategy to be adopted. From a scientific point of view, this should be reflected in the use of methodologies with limited complexity, rapid applicability and low dependence on computing time, as well as in the transparent availability of the methodologies used.

Methods for estimating tolerable by-catch rates should be kept simple, otherwise they may not be understood by all stakeholders and may therefore be difficult to apply.

## The common dolphin, an elusive species

The notion of elusive species covers both demographic and management aspects. From a demographic point of view, these are species for which it is difficult to obtain data and vital rates. It is difficult to get an idea of the viability of the population on the basis of vital rates. It is also a species whose conservation is important in terms of the ecosystem balance (i.e. keystone) and therefore in terms of international legislation, especially in Europe. In Europe, the population is difficult to manage. Limited demographic knowledge is a contributing factor (Murphy et al. 2021). But difficult management discussions between fisheries, politicians, scientists and NGOs are also a contributing factor. Data on bycatch of this species in French waters are also very limited. The European Commission has set up an Infringement Procedure against France in July 2020:

"The Commission is asking France, Spain and Sweden to implement the measures required under the Habitats Directive (Council Directive 92/43/EEC) and Common Fisheries Policy to avoid unsustainable by-catches of dolphin and porpoise species by fishing vessels. [...] Despite well-documented evidence of these species being caught in fishing gear, the problem persists. France, Spain and Sweden have not taken sufficient action to monitor by-catches in their waters and by their fleets, nor made full use of the possibilities that the Common Fisheries Policy provides to comply with their obligation under the Habitats Directive and protect these species. [...] As France, Spain and Sweden have not taken the necessary steps to address these deficiencies, the Commission is sending letters of formal notice to the three countries, which will have three months to address the shortcomings raised. Otherwise, the Commission may decide to send a reasoned opinion"¹.

One of the difficulties in monitoring this population is that it is difficult to define and relate it to the OSPAR Management Unit. The next section is a discussion of the definition of this population.

<footnote>1: Follow this link to see the full procedure of the European Commission: https://ec.europa.eu/commission/presscorner/detail/en/INF_20_1212</footnote>
## 7.2 Definition of the North-East Atlantic common dolphin population

### Population migration and possible distributional changes

#### Exchanges between stocks and management units

Animal populations can exchange individuals through the phenomenon of migration which may occur between different stocks in the case of the common dolphin. Migration is composed of immigration which is the arrival of new individuals, including breeders, into the population, and emigration which is the final departure of individuals from the population. Sometimes the flow of individuals between two populations is not reciprocal. In this case, we speak of a source population that provides individuals to a sink population that receives them. This migration flux implies that the ecological conditions are met. Firstly, there must be a high production of individuals in the source population. Secondly, the sink population must suffer a high loss of individuals (e.g. through high extrinsic mortality) leaving vacant places to occupy the trophic niche of the species. Thirdly, there must be a migration corridor between the two populations with no geographical barrier.

Stocks can exchange individuals through dynamic source (donor) and sink (receiver) migration processes.

<figure><img src="image_125.png" /><figcaption>**Figure 7.4: Possible management units for the North-East Atlantic common dolphin.** The current management unit recognises a panmictic population by taking into account the Bay of Biscay, Iberian coast, Celtic seas and North sea. This management unit could be divided into three separate units. It could also be extended to a larger scale, including much of the North Atlantic.

Figures from Perrin (2018), ICES (2021a) and Caurant et al. (2009).</figcaption></figure>

The vital rates obtained in this project are not consistent with the abundance estimates, which raises the question of the definition of the management unit. The vital rates are associated with a large majority of dolphins presumed dead on the Bay of Biscay continental shelf. Population projections show extinction in 50 years, as do the projections made by Mannocci et al. (2012). The population status also appears to have deteriorated since 1997, with an increasingly low growth rate. But these estimates are not consistent with the observed trend in abundance. The population

There is an inconsistency between the population projections and the observed trend in abundance in the Bay of Biscay and in the Management Unit.
does not appear to be declining at the scale of the Bay of Biscay and the management unit covering the Bay of Biscay, Iberian coast, Celtic seas and North Sea. It is possible that the Bay of Biscay shelf is associated with a sink population that receives individuals from other continental shelves (e.g. Iberia and Portugal) and Celtic seas. In this case, it is possible that the currently recognised management unit is too large and needs to be redrawn according to the classification proposed by Caurant et al. (2009).

Another possibility is that the resulting vital rates are representative of the management unit recognised by OSPAR but that this management unit receives external individual inputs allowing it to be stationary. It is assumed that the population is not viable given current estimates of bycatch numbers across the management unit (ICES 2019). But the population abundance seems to be stationary at the scale of the management unit (Hammond et al. 2017). The currently recognised management unit may not include the majority of the common dolphin population. Migratory flows between the currently recognised management unit and other stocks could explain the stationarity of abundance. This hypothesis is difficult to verify. But it is possible that the currently recognised management unit is not the most appropriate one (Figure 7.4). On the one hand, by omitting external migrations. On the other hand, by not being sufficiently fragmented.

The management unit may not be fragmented enough. Several sub-units should be defined as there may be sink-source exchanges between continental shelf stocks.

The management unit may not include the majority of the common dolphin population. In this case, it may not be large enough and should be expanded to the North Atlantic scale.

<figure><img src="image_126.png" /><figcaption>**Figure 7.5: Common dolphin sightings from the SAMM surveys.** SAMM surveys were conducted in Winter (November 17th 2011 to February 12th 2012) and Summer (May 16th to August 8th 2012). There were more common dolphins observed in summer on the oceanic part of the Bay of Biscay than in winter.

Figure from (Laran et al. 2017)</figcaption></figure>

It is difficult to have a precise idea of the distribution of the common dolphin within the management unit but it is possible that seasonal variations in abundance in the Bay of Biscay are linked to source-sink dynamics. Abundance between summer and winter on the slope and the oceanic part of the Bay of Biscay appears to be different (Figure 7.5). The influx of individuals in summer
may be linked to the breeding season (Murphy et al. 2021, 2009; Murphy, Collet, and Rogan 2005). These seasonal changes have been identified by the Seasonal Abundance of Marine Mammals (SAMM) surveys (Laran et al. 2017). Common dolphin sightings are from winter (17 November 2011 to 12 February 2012) and summer (16 May to 8 August 2012). The winter period sampled mixes months of very low density (November-December) and months of very high density (February) (unpublished data). In my opinion, one main highlight of these surveys is that the oceanic part of the Bay of Biscay may be subject to immigration during Summer (reproduction period). The abundance of small-sized delphinids (common dolphins and striped dolphins) in the Bay of Biscay in Summer 2012 was estimated at 493,591 individuals (95% CI: 342,000 - 719,000) whereas in winter 2011 it was estimated at 284,894 (95% CI: 174,000 - 481,000). It is possible that the difference in abundance of individuals is due to the contribution of non-resident individuals in the Bay of Biscay that come to breed in summer. It is difficult to say whether some of these individuals then remain in the Bay of Biscay to live there. It is also difficult to know the sex and age of the individuals that might move.

## Overall change in distribution

These seasonal movements could be nested within larger spatio-temporal phenomena. Oceanographic processes are known to be structured in space and time. One way to represent this structure is the Stommel diagram generally involves two dimensions (space and time) but may involve a third (e.g. biomass). These representations are initially used to characterise the spatio-temporal phenomena of marine ecosystems (Haury, McGowan, and Wiebe 1978; Margin-Figure 7.6). These diagrams are useful for gaining a better understanding of the spatio-temporal phenomena structuring a system and their scales. They have also been applied to characterise the development of the biosphere (Holling 1992) and also the socio-ecological systems (Holling, Gunderson, and Ludwig 2002). I here provide a rough example of diagram that can be used as an approximate representation (Figure 7.7). The advantage of such a representation is that it allows the linking of different potentially related phenomena that do not take place at the same spatiotemporal scales. Here, it allows us to move from changes in the life history of a resident population, which may be driven by migration (when allocating resources for breeding or feeding), to larger scale phenomena such as distributional changes.

The presence and residence of top predators is strongly influenced by the presence of their prey. Low prey availability can

There seems to be a difference in dolphin abundance between summer and winter in the Bay of Biscay. This difference mainly concerns the continental slope and the oceanic part of the Bay of Biscay.

<figure><img src="image_127.png" /><figcaption>**Figure 7.6: Marine biology Stommel Diagram.** Emphasizes factors in marine biology, with an emphasis on biological productivity, here labeled "biomass variability". Figure from Haury, McGowan, and Wiebe (1978)</figcaption></figure>
<figure><img src="image_128.png" /></figure>

**Figure 7.7: Rough Stommel diagram for common dolphin.** This raw dia-
gram presents processes associated
with the population at three spatio-
temporal scales. The population can
be resident of the Oceanic or Shelf
parts for about two months. Individ-
uals may migrate within the Bay of
Biscay each year. The population may
change its distribution on a larger spa-
tial and temporal scale later. This dia-
gram is only a rough and conceptual
view and should not be considered
as completed.

influence the energy budget of the predator and ultimately influ-
ence its demography (MacLeod et al. 2007). The harbour porpoise
(*Phocoena phocoena*, Linnaeus 1758) is the smallest odontocete in
Europe with a length of about 1.5m and a weight of 50kg; it is also
found in the North-East Atlantic ocean. This species has a shorter
lifespan, a younger age at first reproduction and a shorter calving
interval. Because of its small size, the harbour porpoise cannot
survive long periods without feeding and almost dies after 3 days
without food (Kastelein et al. 1995). The harbour porpoise was
widespread in the North Sea and fed mainly on sandeels (genus
*Hyperoplus*). These fish species are essential to the ecosystems of
the North Sea (Greenstreet 1996). They are also essential for the sur-
vival of harbour porpoises in spring. This season is associated with
the coldest water temperatures in the North Sea, which requires a
high energy intake to maintain a stable internal body temperature
(Watts, Hansen, and Lavigne 1993). The presence of sandeels is
also crucial for the survival of weanlings that are about to take
their first food intakes.

The viability of the predator pop-
ulation depends partly on its prey.
The harbour porpoise, for example, is
heavily dependent on sandeels in the
North Sea.

<figure><img src="image_129.png" /></figure>

**Figure 7.8: Density of Harbour porpoise estimated from the SCANS and CODA campaigns.** The time period under consideration ranges from 1994 to 2016 with 3 main estimates. The density of harbour porpoises in the North-East Atlantic appears to have possibly changed over the period. The population now seems to be denser in the Western Channel, Celtic seas and Bay of Biscay than before.

Figures from Hammond et al. (2017, 2013).