Population assessment of tropical tuna based on their associative behavior around floating objects

Estimating the abundance of pelagic fish species is a challenging task, due to their vast and remote habitat. Despite the development of satellite, archival and acoustic tagging techniques that allow the tracking of marine animals in their natural environments, these technologies have so far been underutilized in developing abundance estimations. We developed a new method for estimating the abundance of tropical tuna that employs these technologies and exploits the aggregative behavior of tuna around floating objects (FADs). We provided estimates of abundance indices based on a simulated set of tagged fish and studied the sensitivity of our method to different association dynamics, FAD numbers, population sizes and heterogeneities of the FAD-array. Taking the case study of yellowfin tuna (Thunnus albacares) acoustically-tagged in Hawaii, we implemented our approach on field data and derived for the first time the ratio between the associated and the total population. With more extensive and long-term monitoring of FAD-associated tunas and good estimates of the numbers of fish at FADs, our method could provide fisheries-independent estimates of populations of tropical tuna. The same approach can be applied to obtain population assessments for any marine and terrestrial species that display associative behavior and from which behavioral data have been acquired using acoustic, archival or satellite tags.


Appendix 2 Details on the models employed to fit the survival curves
The survival curves of residence and absence times were fitted using the three following equations: 46 S(t) = f 1 e −k 1 t + (1 − f 1 )e−k 2 t (S2) and S(t) = β α /(β + t) α The single exponential model in equation ((S1)) corresponds to a time-independent, memoryless process characterized by a single timescale 1/k 1 , where k 1 denotes the constant probability for a failure event to occur. The double exponential model in equation ((S2)) characterizes a time-independent, memoryless process associated to two timescales 1/k 1 and 1/k 2 , with two constant probabilities k 1 and k 2 . Here, f 1 represents the proportion of events associated to the timescale 1/k 1 . Finale the power-law model in equation ((S3)) describes a time-dependent dynamics, with α being the power coefficient and β being the minimal time for the first failure event to occur. In such model the probability for a failure event to occur decays with time as α/(β + t).

Appendix 3 Double exponential model for continuous residence times
Consider a fish population of size N, where the associated fish shows two possible behavioral modes relative to the time spent at the FADs, S and L (short and long residence times). The total number of fish associated with FAD i can be expressed as: where L i (S i ) represents the amount of fish in state L (S) at FAD i. The total population is: where X u represents the number of unassociated fish. If the probabilities to join or depart from a FAD are time-independent, the number of associated fish in each behavioral state at FAD i is described by the following association dynamics: Considering equation ((S6)) at equilibrium leads: From the above equations, the total number of associated fish can be expressed as: S8) and the ratio between the associated and total number of fish is given by: The parameters on the r.h.s. of equation (S9) can be inferred from the survival curves of CRTs and CATs. In the presence of two behavioral states, the CRTs recorded at FAD i follow a double exponential model: Where C S i represents the proportion of residence times for the behavioral state S and θ S i and θ L i are the two probabilities to depart from FAD i for a fish in behavioral state S and L, respectively. Reversely, the survival curves of CATs follow a single exponential model of the form: (S10), the values of θ S i and θ L i can directly be inferred from the fit of the survival curves of CRTs. Moreover, the probabilities µ S i and µ L i to reach FAD i for each behavioral state can be expressed as: where n i is the number of CRTs recorded at FAD i and n tot = ∑ i n i and C S i and µ tot can be inferred from equations (S10) and (S11), respectively. Finally, considering the limit θ S i θ L i , equation (S9) can be reduced to a single exponential model of equation (5), where the only relevant timescales are related to the long residence times θ i = θ L i and µ i = µ L i .

Survival curves of CRTs for FAD-class 1
The lowest values of the AIC were found for the double exponential model, see Table S5. However, one of the exponents of the double exponential function (k 1 ) was not significantly different from zero (p-value=0.16). As a consequence, the double exponential model was rejected. The AIC of the single exponential and the power-law models were very close, but the standard errors of the power law model were very high. For model parsimony, we therefore considered the single exponential model as the best fitting function for FAD-class 1.

Survival curves of CRTs for FAD-class 2
The double-exponential model showed the lowest AIC, which was significantly lower than the other models (Table S5). For this reason, this model was considered in the rest of the analysis.

Survival curves of CATs
The AIC of the single and double exponential models were very close (Table S5). However, one of the exponents of the double exponential function (k 2 ) was not significantly different from zero (p-value=0.7). As a consequence, the double exponential model was rejected. Reversely, the power-law model could not be fitted to the data, despite multiple trials were conducted by considering different initial conditions.

Appendix 5 Details on the model application to Yellowfin tuna
From the results obtained from the survival analysis of CRTs and CATs (see Results section and Appendix 4), the estimate of the abundance ratio took into account: • An heterogeneous FAD array at equilibrium, with two classes of FADs, the fist class (FAD-class 1, with one FAD only) being characterized by a single exponential model and the second class (FAD-class 2, with 12 FADs) by a double exponential model.
• A single timescale related to a single exponential model for CATs.
In this case, the ratio between the associated and total number of fish can be analytically derived by combining equation (5) for a single exponential model and (S13) for a double exponential model for CRTs, leading to: with: where µ 1 , µ S 2 and µ L 2 (θ 1 , θ S 2 and θ L 2 ) represent the probabilities to associate with (depart from) the FADs of Class 1 and Class 2, respectively, considering for the latter two behavioral modes (S and L), see Appendix 3. The departure probabilities θ 1 , θ S 2 and θ L 2 could be estimated from the fits of the survival curves of CRTs (see equation (8) and ((S9))). Similarly, the arrival probabilities µ 1 , µ S 2 and µ L 2 could be estimated from field data through the following equations: where n 1 (n 2 ) is the number of CRTs recorded at FAD-class 1 (2), n tot = n 1 + n 2 is the total number of CRTs, C S 2 is the fraction of CRTs associated to short residence times for FAD-class 2 (equation (S8)) and µ tot is the probability to associate with one of the FADs of the array related to the survival curves of CAT (equation (S10)). In the case where CRT1 was not considered, the number of fish tagged at the FAD of each class (Table S1) were subtracted from the estimated values of n 1 , n 2 and n tot and the probabilities of joining each FAD class recalculated from equation (S15). Table S6 resumes the values of n 1 , n 2 and n tot with and without considering CRT1. Similarly, the approximate formula for the abundance ratio in the limit θ L 2 θ S 2 could be written as: The model parameters used in the stochastic simulations to reproduce the observed association dynamics are shown in Table S7. The model parameters of FAD class A correspond to those of FAD-class 1 (FAD HH) in the field data, whereas FAD-class B denotes the other FADs. The probabilities to reach each of the FADs of FAD-class B were derived from equation (S15) by dividing the estimated values of µ L 2 by the number of FADs that constituted that class (12 FADs). Table S1. Tagging strategy: FAD of tagging, year and month of tagging, number of yellowfin tuna tagged and size range (mean fork length ± SD).
Feb Mar Apr May June July

5/6
Parameter symbol -Name Value Value (*) N -Total number of fish 1.0e4 1.0e4 N T -Number of tagged fish 10 10 p A -Total number of FADs in class A 1 1 p B -Total number of FADs in class B 12 12 µ A -Probability to reach FAD-class A 0.0850 0649 µ B -Probability to reach FAD-class B 0.0174 0.0184 θ A -Probability to depart from FAD-class A 0.047 0.047 θ B -Probability to depart from FAD-class B 0.27 0.27 T start -Time of tagging 1.0e4 1.0e4 T end -End Time of simulation 1.0e5 1.0e5 Table S7. Model parameters for reproducing the experimental data. The shaded cells represent the model parameters that are obtained from the experimental data (Table 3). Last columns denoted with (*) correspond to model parameters with probabilities of joining the FADs obtained when excluding CRT1. Notice that the probability to reach a single FAD of FAD-class B is obtained by considering µ L 2 in equation (S15) and dividing by the total number of FADs of class B p B . The values of the other parameters are the same as in the stochastic simulation.