Deterministic logistic model
The following population dynamics model was applied to reconstruct the initial dugong population size in 1894 from fishery statistics between 1894 and 1914:
$$N_{t + 1} = N_{t} left( {1 , + r{-}r , N_{t} /K} right) – C_{t} ,$$
where r is the intrinsic rate of population increase, Nt is the population size in year t, K is the carrying capacity, and Ct is the number of individuals removed from the waters near the Ryukyu Islands in year t. The carrying capacity (K) in 1893 was sufficient to sustain the initial population of dugongs at that time (N1894). The intrinsic rate of population increase (r) was given between 1 and 5% within a range of natural one.
Approximate Bayesian calculation
We conducted approximate Bayesian calculation (ABC)32 to estimate the number of individuals in 1979 based on bycatch data between 1979 and 2019, and the constraints of the numbers of individuals were 11 in 1997, three in 2007, and almost extinct in 2019. We denoted fecundity as f, the survival rate until 1 year old as s0, the annual survival rate after 1 year old as s, the age at maturity as am, and the physiological longevity as A. We assumed that the sex ratio at birth was 1:1 on average; the age at maturity am was eight years of age33, and the physiological longevity A was 73 years6. We ignored environmental stochasticity because no mass deaths caused by infectious diseases or changes in survival or mortality rates due to environmental fluctuations have not been recorded during this period. We also ignored density effects because the carrying capacity of the location was sufficiently greater than the initial population size, and our goal was to investigate the possibility of population recovery after a decrease in population using a population dynamics model and estimate the natural growth rate during this period. The detailed extinction risk depends on age structure.
According to the life history parameters, except the physiological longevity compiled by (ref.33), the annual survival probability of an a year-old individual is s for a = 1, 2, …, 72; s0 for a = 0, and 0 for a = 73; the reproductive probability of an adult female > 8 years old is 2f. As the number of years for a population to become extinct or recover depends on age composition, age-specific survival, and reproductive rates, we obtain the population growth rate by the maximum eigenvalue of the following Leslie matrix, L = {Lij} (i = 1,…73, j = 1,…,73) as:
$$L_{i1} = s_{0} f/2quad {text{for}}quad i ge a_{m} ,L_{i+ 1,i} = squad {text{for}}quad i = 1, ldots ,72,quad {text{and}}quad L_{ij} = 0,{text{otherwise}}{.}$$
We used the population growth rate λ, defined by the maximum eigenvalue of L, as an indicator of the population growth rate.
We assumed that the sex of each individual in 1979 was randomly sampled by the 1:1 sex ratio, and its age was randomly sampled by the stable age structure that is given by the eigenvector of the Leslie matrix with the maximum eigenvalue. We assumed that the number of individuals at age 1 year in year t + 1, denoted by N1,t+1, is determined by the binomial distribution:
$$Prleft[ {N_{1,t + 1} = x} right] = left( {begin{array}{*{20}c} {N_{f} } x end{array} } right)left( {s_{0} f} right)^{x} left[ {1 – left( {s_{0} f} right)} right]^{{N_{f} – x}} ,$$
where Nf represents the number of adult females in year t. We assumed that no twins were born. We assumed that the probability that an individual with age x survived in the next year is s if x = 1 or s0 if x = 0. We also assumed that Ct individuals who died by bycatch were randomly chosen from any sex and age because the age of individuals caught by bycatch is rarely known. We do not know the sex of some individuals.
We assumed the following prior distributions for N1997, f, and s: N1979 (in) U(11, 80), f (in) U(1/14, 1/6) if at least one adult male existed in the population, s0 (in) U(0.1, 0.85); and s (in) U(0.8, 0.97), where U(a, b) is the uniform random variable between a and b. These probabilities were constant for each simulation trial from 1997 to 2019. We selected the set of parameters with the population growth rate (λ) obtained when the maximum eigenvalue of the Leslie matrix was between 0.96 and 1.01.
We rejected trials that did not satisfy the following summary statistics: N1997 ≥ 11 (intensive survey in 1997), Nt ≥ 3 during 2004–2017 (monitoring), and N2019 ≤ 1 (“local extinction”). We obtained the prior distributions of N1997, f, s0, s, and N2004, and of the > 130,000 trials in the prior distribution with natural population growth rates λ of 96.1–98.8%, 99.3% were rejected. For 95% of the 1000 adopted trials, N1979 ranged from 14 to 58. If λ > 98%, N1997 was ≤ 45 for the adopted trials (Extended Data Fig. 7. Even if all the stranding deaths were due to anthropogenic factors, such as the release of dugongs after bycatch or boat strike, the range of N1997 changed to < 56 if λ > 98%, with only a slight upward shift, but positive natural growth rate (or λ > 1) was again very unlikely (0.3%) among the adopted trials.
Population viability analysis to assess the impact of bycatch on the extinction risk
We re-evaluated the extinction risk with and without bycatch using the 1000 parameter sets of N1979, f, s0, and s that satisfied the summary statistics in the ABC and stochastic individual-based model, beginning from N1979 for the corresponding parameters. For each parameter set, 100 trials were conducted for each scenario to compare the extinction risks.
Source: Ecology - nature.com
