Incorporating Stakeholder Knowledge into a Complex Stock Assessment Model: The Case of Eel Recruitment

: Mistrust between scientists and non-scientist stakeholders is a key challenge in ﬁshery management. The support of management with complex models is difﬁcult, as these models cannot easily communicate their results to end users, resulting in a lack of conﬁdence. As an example, the complex life cycle of the European eel raises problems of coordination and discussion among the multiple actors involved in the species’ management. The GEREM model has been proposed as a tool for estimating recruitment, but its complexity, which is essential for addressing the characteristics of the species, makes it difﬁcult to understand and accept by all stakeholders. In the context of the SUDOANG project, we proposed a co-parametrisation of this assessment model to tackle this mistrust. Through the use of various questionnaires, stakeholders were involved in two important choices for the model (zone deﬁnition and prior construction). Regular workshops and presentations were organised to explain the model rationales and to gather feedback and expectations. The results show that stakeholders have very similar perceptions of the potential deﬁnitions of sub-areas of recruitment in south-western Europe, and these perceptions are consistent with the underlying environmental conditions. On the other hand, the stakeholders have contrasting opinions about the exploitation rates of ﬁsheries in different river basins, and the use of their knowledge currently has little effect on GEREM estimates. More importantly, the overall approach of this study is thought to have reconstructed the trust and conﬁdence among participants.


Introduction
There is a long history of misunderstandings between scientists and non-scientific actors in fisheries management [1], and rebuilding trust among them is a key challenge [2][3][4]. Fisheries dynamics are hardly predictable [5] and display chaotic behaviour [1,6]. While scientists have developed increasingly complex numerical models and seek general rules three spatial scales: the south-western Europe scale, sub-region scale and river basin scale. By doing so, GEREM provides insight into the overall trend in recruitment (on the largest spatial scale). By providing estimates on a fine spatial scale (e.g., river basins), GEREM also allows local managers to compare recruitment estimates with escapement estimates provided by EDA and, consequently, to assess trends in mortality rates during the growth phase of eels in the water bodies under their jurisdictions. The use of these two models is a step forward compared to the current situation in which managers of the Sudoe area use very different data and methods to respond to the Eel Regulation, weakening the consistency of the estimates and impairing comparisons of the results. However, a certain level of complexity is required in the two models to address the complexity of the life cycle of the species. Therefore, to achieve the SUDOANG objectives, there was a need to build stakeholder confidence in the models.
In this article, we present how the GEREM model was implemented during the SUDOANG project and, more specifically, how knowledge of stakeholders was elicited and used in the model tuning, with the underlying assumption that explicitly using their knowledge was a way to acknowledge their expertise and to increase their interest and confidence in the stock assessment model.

GEREM Model: Assumptions and Data Requirements
The GEREM model is described in detail by [21,28]. Since the aim of this paper is not to discuss the results of the model, but to focus on the use of expert knowledge and on its influence on the result, we just recall the main underlying assumptions and describe the modifications made for this exercise. GEREM is a Bayesian state-space model that estimates yearly recruitment on three nested spatial scales: river basin scale R c,z (y), sub-region (zone scale) R z (y) and overall study area.
Mean recruitment in a river basin is assumed to be proportional to the corresponding zone recruitment with the proportionality factor equal to a function of the surface of the river basin. This assumption allows using recruitment estimates in some basins to extrapolate recruitment to other basins of the same zone. However, for this extrapolation to be reasonable, the study area must first be divided into regions, groups of river basins, where recruitment is thought to be homogeneous in terms of temporal trends and recruitment levels.
To estimate recruitment, GEREM is fitted on the available time series of recruitment that are collected in various river basins. The model assumes that an observed time series of abundance I A i,c (y) collected in a river catchment c follows a lognormal distribution with a mean proportional to the model-predicted abundance: Most time series are relative time series of abundance; i.e., they only provide information on the trend of recruitment but not on the actual abundance because no information is available on the scaling factor q i . GEREM fitting requires some time series that provide information on the real recruitment in some basins (at least one series per zone) to extrapolate recruitment in any river basins. Different kinds of time series can be used [29]. In some river basins, models have been used to estimate absolute recruitment (e.g., [30][31][32]). In such time series, the scaling factor q i is thought to be equal to 1. If the exploitation rate is constant, commercial catches can also be used as time series of abundance. In that case, the scaling factor is necessarily less than 1 (as catches are necessarily less than recruitment), and more precise information can be used through informative priors. The same line of reasoning exists for time series of glass eel ladder counts; scaling factors that correspond to the overall trap efficiency are necessarily less than 1. Time series for which information on the scaling factor q i is available (commercial catches, absolute estimates and trap monitoring) provide insight into the absolute level of recruitment R c,z (y) which can then be extrapolated to other river basins in the same zone and by doing so, information on zone recruitment R z (y) and on overall recruitment R(y) can be obtained.

Stakeholder Involvement in the SUDOANG Project
The programme Sudoe V, and therefore the SUDOANG project, covers the Iberian Peninsula and two French regions (Nouvelle-Aquitaine and Occitanie). The GEREM modelling exercise was slightly extended and ranged from western Brittany to south-eastern France. The project gathers a large panel of stakeholders from three countries: scientists, commercial fishers, administrators in charge of species protection and/or habitat protection, and control forces in charge of the fight against the illegal trade of protected species (see Table S1 in supplementary material). To promote the co-construction of tools, different meetings, workshops and online surveys were organised during the project to keep partners informed about the work progress, to obtain their feedback and to determine their expectations. More specifically, regarding GEREM, the following programmes were organised:  Figure S1). • A workshop was then organised (15-16 November 2018) to discuss and validate the recruitment zone definition and plan an exercise to collect information from the participants regarding the exploitation rates of commercial fisheries. The participants were each given a sheet with the main glass eel fisheries and asked to provide a value and a range. • A document summarising the main findings about the work on exploitation rates was sent a few months later, and an online questionnaire was sent to all participants of the project to re-run the exercise.

•
The results regarding the opinions on exploitation rates were presented in the second annual meeting (19- For all questionnaires, participants were allowed to remain anonymous so that they felt completely free to express their opinions.

Available Data Series
As with all models, a data inventory was required before fitting the model. GEREM can use any kind of time series of recruitment as long as the data are considered to be proportional to the real abundance, i.e., a reliable abundance index. Moreover, GEREM needs time series for which information (even vague information) on the scaling factor q i (Equation (1)) is known; therefore, it is important to list the type of time series used. An inventory of available time series was made during the kick-off meeting in Lisbon (Table 1 and Figure 1). GEREM aims to estimate recruitment (since European eel is a catadromous species recruitment, recruitment corresponds to the quantity of glass eels that enter river basins), selected time series were collected as downstream as possible, before being affected by any local anthropogenic pressures (e.g., fishery, obstacles). Most of the time series had already been used by [21] and are also used by the WGEEL [23], and descriptions of the time series are provided in these references. The only new series was a commercial catch time series, which corresponded to the daily catch from an artisanal commercial fishery in the Ter River (eastern Spain). While this time series is short (only 4 years), it was possible to provide information on the corresponding exploitation rate (q i ). A description of the fishery is available in the Catalunyan Eel Management Plan [33]. Table 1. Available time series of recruitment and key characteristics. Type indicates whether a time series is purely relative or if it is possible to provide information on the scaling factor (Equation (1)). Surface corresponds to the surface area of the river basin. The length of the time series is indicated by its first year, its last year and the number of data points (Nb data). Their locations are visible on Figure 1. 1 and Figure 1). GEREM aims to estimate recruitment (since European eel is a catadromous species recruitment, recruitment corresponds to the quantity of glass eels that enter river basins), selected time series were collected as downstream as possible, before being affected by any local anthropogenic pressures (e.g., fishery, obstacles). Most of the time series had already been used by [21] and are also used by the WGEEL [23], and descriptions of the time series are provided in these references. The only new series was a commercial catch time series, which corresponded to the daily catch from an artisanal commercial fishery in the Ter River (eastern Spain). While this time series is short (only 4 years), it was possible to provide information on the corresponding exploitation rate (qi). A description of the fishery is available in the Catalunyan Eel Management Plan [33].

Definition of Model Zones Using Questionnaires and Consistency with Environmental Conditions
After listing available time series, a second task before fitting the model is to define the zones. The model assumes that yearly recruitment in basins within a given zone follow similar temporal trends and that recruitment from one basin can be extrapolated to another basin within the same. Finally, we needed at least one time series for which information on the scaling factor is available. In [28], the authors used French EMUs to define regions, whereas [21] used scientifically defined ICES ecoregions. Here, we decided to gather stakeholder opinions about the zone definitions. After the presentation of the model and available time series and after having clearly explained the assumptions that should be met, participants of the kick-off meeting were asked to draw on a map the zones they would use given their knowledge. A map corresponded to the clustering of all the rivers basins of the area (i.e., not only rivers basins in which time series are available) into a variable number of zones ( Figure S1). The number of zones was unrestricted. Based on each map, it was possible to determine whether two river basins were classified in the same zone according to expert opinion.
The DISTATIS method [34,35], applied with the DistatisR package [36], was used to produce maps where the distance between points reflects their similarity as in PCA (principal component analysis, with each point standing for a river basin. It includes the following steps: (1) for each map drawn by a stakeholder, a similarity matrix was constructed, with one row and one column per basin, that holds 1 if two basins were classified as being within the same zone and set to 0 otherwise. (2) Then, theses matrices were aggregated into an optimal compromise matrix, which is a cross-product matrix (similar to the correlation matrix in a PCA) quantifying how distant pairs of river basin, and represents an optimal consensus among stakeholders. (3) This consensus matrix can then be analysed using an eigen decomposition (just as in a PCA), which produces a "compromise space" of river basins. (4) Then, individual expert classification matrices can also be projected onto this consensus space to analyse similarities among experts. Two types of diagrams can be produced that are similar to PCA: • A plot of experts in which each point stands for an expert and the distance between points stands for the distance between the opinions of the experts. This diagram is not produced here, as many stakeholders remained anonymous. • A plot of river basins on which each point stands for a river basin and two points are close when most experts thought they should be in the same zone. For simplicity, only the first two components were plotted. However, hierarchical clustering was carried out to classify river basins based on their distance using all the coordinates of all dimensions from the PCA.
To check whether the resulting zones were consistent with environmental conditions, a second type of analysis was independently carried using only environmental conditions instead of expert knowledge ( Figures S2 and S3). This analysis was based on the following environmental data collected for each river basin:

•
Marine conditions (averaged from 2006 to 2016 for the period from October to April, corresponding to the main period of recruitment [38]), including the average sea surface temperatures (source: IFREMER), average salinity (source: IBI MFC), average concentration of Chlorophyl a (source: ACRI-ST) within a 100 km buffer around the outlet, and the distance to the 150 m isobath as a proxy of the local continental shelf width.
Using these environmental variables, we carried out a through a multivariate regression tree [39]. The environmental variables of each river basin (previously mean deleted and scaled by the standard deviation) were used as explanatory variables of a continuous variable defining their positions: the coastline distance of each river basin outlet to Monaco, Water 2021, 13, 1136 7 of 24 located at the extreme south-east of the study area. The tree was fitted using the partykit package [40]. Regression trees seek successive criteria/thresholds on the explanatory variable that minimise the dissimilarity within the resulting clusters.
The resulting classification, based on physical and environmental variables were finally compared to the maps built on expert knowledge. This was a way to ensure that the two methods run independently (expert-based classification versus environmentally based classification) and provide consistent results.

Construction of Priors on Glass Eel Fishery Exploitation Rates
As explained in Section 2.1, knowledge on the scaling factors of time series is valuable if it is available. Among the available time series, 8 time series corresponded to commercial catches; therefore, it might have been possible to build Bayesian priors on the scaling factors of these time series. This would ensure that the model explored zones that are thought to be reliable by stakeholders.

Design of the Questionnaires
This task was developed in two stages. During the workshop in November 2018, after a presentation detailing the progress of the work, we introduced the notion of exploitation rates and how the knowledge of participants could be used in the model. We presented this line of reasoning: "if the model estimates that the fishery harvests 100% of incoming glass eels in a given estuary, you would say it is impossible. Conversely, if it estimates an exploitation rate of 0%, you would also say that this is impossible. We would like to know the minimum and maximum bounds that you consider to be possible in this estuary". The idea here was not to ask for precise estimates of exploitation rates but rather to eliminate impossible value ranges, i.e., estimates that would be considered impossible by the stakeholders and would consequently reduce their confidence in the model. To do this, we built a questionnaire in which, for each time series, we presented some characteristics of the considered estuary (especially its width) and of the fishery (the number of gears, size of the gears, duration of the legal fishing season) and asked attendees to provide a minimum bound and a maximum bound of the plausible exploitation rates in each estuary with the clear instruction of focusing on the elimination of impossible values, not on giving a precise estimate. After this face-to-face workshop, a report was made to summarise the approach and the main results of this first trial to both the attendees of this workshop and to other SUDOANG participants who had not attended the meeting. A similar questionnaire was then proposed online to facilitate the participation of as many partners as possible (see an example in Supplementary material).

Types of Priors
To see the effect of incorporating expert knowledge on the model outputs, we constructed three types of priors for the exploitation rates of commercial catch time series, namely, "control", "uniform" and "triangular". As a "control experiment", we first used totally uninformative priors that did not consider stakeholder opinions at all, using uniform distributions between 0 and 1; these were called the "control" priors: For the two other types of priors, which were based on stakeholders' opinions, we followed the approach of [41]. For each time series of commercial catch i and each stakeholder having expressed an opinion, we drew 1000 values of the scaling factor in a uniform distribution ranging from the minimal bound to the maximum bound given by the stakeholder. Then, for each time series, a beta distribution was fitted on the mixture of all drawn samples. This resulted in a so-called "uniform" prior of the following form: 8 of 24 where MI N i (MAX i ) is the overall minimum (maximum) bound provided by stakeholders for the corresponding exploitation rate. We also tested an alternative approach to build a more informative prior. The approach was very similar to the previous one, except that samples for each time series and stakeholder were drawn from a triangular distribution instead of a uniform distribution. A triangular distribution is a triangular-shaped density distribution in which the density is zero outside the basis of the triangle. For each expert, we used their reported minimum and maximum bounds to define the basis of the triangle, and the peak of the triangle was defined as the average value between the minimum and maximum bounds. Then, similar to the previous approach, we drew 1000 scaling factors within this distribution and fit a beta distribution on the mixture of all drawn samples. This approach was used to build "triangular" priors:

Running the Model
The other model settings were identical to those used by [21], except for the scaling factor for the prior time series for which we used a more vague prior: which corresponded to a vague prior between 0 and 1/3 (eel fishway passabilities are generally estimated to be approximately 1/3) [42,43]. However, all eels do not reach the fishway, so the count multiplied by 3 should be considered a minimum value for recruitment in the basin.

Model Fitting
The descriptions of the river basins were based on three national hydrographics networks that were pooled together during the project [44]. For example, the French component was based on the RHT network [45] that had been used by [28] to fit GEREM in France. The hydrographic networks provided the shape and the surface of the basins (parameters that are required in Equation (1)). The study region had an overall surface area of 953,268 km 2 and was made up of 2582 river basins.
The model was fitted with JAGS [46] using the package runjags [47]. For the three types of priors (e.g., "control", "uniform" and "triangular"), 3 chains were run independently with an adaptive phase of 80,000 iterations and a burnin period of 100,000. Then, the chains were run for 50,000 iterations, with a thinning rate of 50 iterations, leading to 1000 iterations kept for inference per chain.

Zone Definition
During the kick-off meeting, 24 maps were produced by the 47 stakeholders attending the meeting (23 scientists, 5 fishermen, 14 managers, 1 NGO and 4 control forces), each map corresponding to a suggested zone definition. Stakeholders proposed from two to nine zones. The proposal with two zones cut the region between the Atlantic and the Mediterranean parts of the study area, a separation upon which all participants agreed.
The PCA (Figure 2) showed three main clusters of river basins. The first axis separated the Mediterranean basins (positive values) from the Atlantic river basins (negative values). The second axis corresponded to a latitudinal gradient on the Atlantic coast. These two first axes represented 59% of the total inertia. nine zones. The proposal with two zones cut the region between the Atlantic and the Mediterranean parts of the study area, a separation upon which all participants agreed.
The PCA (Figure 2) showed three main clusters of river basins. The first axis separated the Mediterranean basins (positive values) from the Atlantic river basins (negative values). The second axis corresponded to a latitudinal gradient on the Atlantic coast. These two first axes represented 59% of the total inertia. The diagram of the inertia gain of the clustering suggested the partitioning of the basins into three to five groups (Figures 3 and S4). Figure 1 shows the map resulting from the partitioning of four groups: the first zone corresponds to the French coast surrounding the Bay of Biscay (ATL_F) and is separated from the second zone corresponding to the Cantabrian Sea (CANT) at Cap Breton Canyon. The third group corresponds to the Atlantic coast of the Iberian Peninsula (ATL_IB), and the fourth group corresponds to the Mediterranean area (MED). The partitioning into three groups consists of merging the CANT and ATL_F zones (leading to a large Bay of Biscay zone). In the partitioning into five groups, the ATL_IB zone is split into two zones at Lisbon. In this last situation, only one relative time series of abundance (GuadG, Figure 1) would have been available for the zone corresponding to the south-western part of the Iberian Peninsula, so this configuration was rejected because the lack of data would not allow for the application of GEREM. The same applied to any number of zones larger than 5.
The maps of environmental variables are presented in the supplementary material ( Figures S2 and S3). The regression tree of these environmental variables led to slight differences from the stakeholder classifications ( Figure 4). The ATL_F zone is almost similar except that its southern border is slightly more to the north in the regression tree than in the expert opinion map. The same stands for the MED zone, whose border lies slightly eastward compared to that in the stakeholder opinion. However, stakeholders were aware of the presence of the Gibraltar Strait and Cap Breton Canyon, two features that were not taken into account in the regression tree of environmental variables. Then, the regression tree produced two additional zones: a zone that corresponds to the CANT zone and integrates part of the ATL_IB zone (the separation between CANT and ATL_IB was retrieved if we extended the depth of the regression tree). Then, the last zone is similar to the 5th zone that was proposed by stakeholders but was not kept because of the absence of available time series. Figure 1 shows the map resulting from the partitioning of four groups: the firs corresponds to the French coast surrounding the Bay of Biscay (ATL_F) and is sepa from the second zone corresponding to the Cantabrian Sea (CANT) at Cap Breton Ca The third group corresponds to the Atlantic coast of the Iberian Peninsula (ATL_IB the fourth group corresponds to the Mediterranean area (MED). The partitionin three groups consists of merging the CANT and ATL_F zones (leading to a large B Biscay zone). In the partitioning into five groups, the ATL_IB zone is split into two at Lisbon. In this last situation, only one relative time series of abundance (GuadG, F 1) would have been available for the zone corresponding to the south-western part Iberian Peninsula, so this configuration was rejected because the lack of data wou allow for the application of GEREM. The same applied to any number of zones large 5.  A dynamic factor analysis [48] confirmed that all available recruitment time series displayed very similar temporal trends (not presented here; [49]), so the assumption of similar trends within zones is not challenged. Therefore, after stakeholder validation, we decided to use the configuration with four zones: ATL_F, CANT, ATL_IB and MED (Figure 1). tree produced two additional zones: a zone that corresponds to the CANT zone an grates part of the ATL_IB zone (the separation between CANT and ATL_IB was ret if we extended the depth of the regression tree). Then, the last zone is similar to t zone that was proposed by stakeholders but was not kept because of the absence of able time series. A dynamic factor analysis [48] confirmed that all available recruitment time seri played very similar temporal trends (not presented here; [49]), so the assumption of s trends within zones is not challenged. Therefore, after stakeholder validation, we d to use the configuration with four zones: ATL_F, CANT, ATL_IB and MED ( Figure 1

Opinions on the Exploitation Rates and Resulting Priors
After the second stage of the questionnaires, 20 stakeholders expressed op about exploitation rates in at least one basin (26 persons attended the presential wor in November 2018 for the first round of questionnaire: 13 scientists, 4 fishermen sentatives, 5 managers and 1 NGO). Their opinions were contrasting ( Figure 5 bounds were large for all rivers, but this was consistent with the objective of the ex of only limiting impossible values rather than providing precise estimates.

Opinions on the Exploitation Rates and Resulting Priors
After the second stage of the questionnaires, 20 stakeholders expressed opinions about exploitation rates in at least one basin (26 persons attended the presential workshop in November 2018 for the first round of questionnaire: 13 scientists, 4 fishermen representatives, 5 managers and 1 NGO). Their opinions were contrasting ( Figure 5). The bounds were large for all rivers, but this was consistent with the objective of the exercise of only limiting impossible values rather than providing precise estimates.
Despite the large bounds in the stakeholder opinions, some differences were visible among the compromise priors ( Figure 6). The Gironde priors were the widest, while the Ter River had priors favouring the lowest exploitation rates. The difference between the "triangular" and "uniform" priors was limited, although "triangular" priors were slightly more informative than "uniform" priors, as seen by the higher probability densities.
Temporal trends in recruitment, model fits to observations and the effect of expertbased priors.
GEREM confirms the decline in recruitment since the early 1980s, regardless of whether stakeholders' opinions are taken into account (Figure 7). All priors provide similar results. However, the stakeholder-based priors lead to slightly higher recruitment estimates: stakeholders-based priors favoured smaller exploitation rates than uniform priors (Figure 6), and since recruitment is assumed to be proportional to the commercial catches divided by the exploitation rate, a smaller exploitation mathematically leads to a higher recruitment.

Figure 5.
Opinions expressed by the different stakeholders regarding the possible values of exploitation rates (many stakeholder answers were anonymous, so they only had an ID). Each segment represents the range defined as "possible" for the corresponding river basin by one stakeholder. An example of questionnaire used to collect the information is presented in supplementary material.
Despite the large bounds in the stakeholder opinions, some differences were visible among the compromise priors ( Figure 6). The Gironde priors were the widest, while the Ter River had priors favouring the lowest exploitation rates. The difference between the "triangular" and "uniform" priors was limited, although "triangular" priors were slightly more informative than "uniform" priors, as seen by the higher probability densities.

Figure 5.
Opinions expressed by the different stakeholders regarding the possible values of exploitation rates (many stakeholder answers were anonymous, so they only had an ID). Each segment represents the range defined as "possible" for the corresponding river basin by one stakeholder. An example of questionnaire used to collect the information is presented in supplementary material.
When analysing the results on the zone scale, several things can be noted ( Figure 8). First, while all zones displayed a decreasing trend in recruitment from the 1980s to the early 2010s, the decrease in ATL_F was maintained to date, while recruitment appeared to stabilise in the other zones after 2010. However, the trend in ATL_F is based on a single time series of recruitment (GirScG) since the upkeep of French fishery-based time series ended after the implementation of the French Eel Management Plan [50]. Similarly, a decrease was estimated in MED during the 1960s; however, these estimates were based only on two fishery-based time series (Ebro and Albufera de Valencia). As before, the use of different priors did not lead to major differences in the temporal trends. However, in the case of ATL_F and ATL_IB, the estimations of recruitment were higher when stakeholder opinions were taken into account. In addition, in the case of ATL_F, the recruitment estimates were higher when using the "triangular" priors than when using the "uniform" priors. Temporal trends in recruitment, model fits to observations and the effect of expert based priors.
GEREM confirms the decline in recruitment since the early 1980s, regardless o whether stakeholders' opinions are taken into account (Figure 7). All priors provide sim ilar results. However, the stakeholder-based priors lead to slightly higher recruitment es timates: stakeholders-based priors favoured smaller exploitation rates than uniform pri ors (Figure 6), and since recruitment is assumed to be proportional to the commercia catches divided by the exploitation rate, a smaller exploitation mathematically leads to higher recruitment. When analysing the results on the zone scale, several things can be noted (Figure 8 First, while all zones displayed a decreasing trend in recruitment from the 1980s to th early 2010s, the decrease in ATL_F was maintained to date, while recruitment appeare to stabilise in the other zones after 2010. However, the trend in ATL_F is based on a singl time series of recruitment (GirScG) since the upkeep of French fishery-based time serie ended after the implementation of the French Eel Management Plan [50]. Similarly, a de crease was estimated in MED during the 1960s; however, these estimates were based onl on two fishery-based time series (Ebro and Albufera de Valencia). As before, the use of di ferent priors did not lead to major differences in the temporal trends. However, in the cas of ATL_F and ATL_IB, the estimations of recruitment were higher when stakeholder opin ions were taken into account. In addition, in the case of ATL_F, the recruitment estimate were higher when using the "triangular" priors than when using the "uniform" priors. The results can also be displayed on the river basin scale and compared with th observations: the fits are very good and are similar among the types of priors (Figure 9).
The posterior distributions of scaling factors qi (Equation (1)), i.e., the exploitatio rates (see Section 2.5), are also similar regardless of the type of priors used, especially fo the Gironde estuary and Adour River ( Figure 10). However, some slight differences wer observed if we compared the posterior distribution obtained with either the "control prior or with the stakeholder-based priors. For example, the probability densities of th very low and high exploitation rates in the Minho basin were almost zero when the stake holder-based priors, i.e., "uniform" and "triangular" priors, were used. In the Albufer de Valencia lagoon, using stakeholder-based priors led to the estimation of lower explo tation rates. The results can also be displayed on the river basin scale and compared with the observations: the fits are very good and are similar among the types of priors (Figure 9).  The posterior distributions of scaling factors q i (Equation (1)), i.e., the exploitation rates (see Section 2.5), are also similar regardless of the type of priors used, especially for the Gironde estuary and Adour River ( Figure 10). However, some slight differences were observed if we compared the posterior distribution obtained with either the "control" prior or with the stakeholder-based priors. For example, the probability densities of the very low and high exploitation rates in the Minho basin were almost zero when the stakeholderbased priors, i.e., "uniform" and "triangular" priors, were used. In the Albufera de Valencia lagoon, using stakeholder-based priors led to the estimation of lower exploitation rates. Figure 10. Posterior density of exploitation rates depending on the type of prior used. The exploitation rates correspond to the scaling factor in Equation (1). Priors were based on expert knowledge (Section 2.5). The control priors are totally uninformative and do not use stakeholder opinions. On the other hand, the two other types of priors are based on stakeholder opinions, among which the triangular priors aim at providing more weights to stakeholder opinions. As in any Bayesian model, the posterior distribution is the distribution of exploitation rates in the MCMC chains of the GEREM model.

Involving Stakeholders in Model Construction
Participatory modelling has been proposed as an important step in overcoming misunderstandings among stakeholders, rebuilding trust and achieving good governance in fisheries management [7]. In this context, this article presents an innovative exercise to involve stakeholders in the parametrisation of a complex model to increase their commitment and confidence in the approach and to achieve the production of a common assessment tool for the European eel in south-western Europe. This question of confidence in complex models as tools used to support fisheries management is crucial, and there is a need to improve transparency regarding model assumptions and uncertainties as well as a need to improve the transfer of model results [12]. The approach used here addressed

Involving Stakeholders in Model Construction
Participatory modelling has been proposed as an important step in overcoming misunderstandings among stakeholders, rebuilding trust and achieving good governance in fisheries management [7]. In this context, this article presents an innovative exercise to involve stakeholders in the parametrisation of a complex model to increase their commitment and confidence in the approach and to achieve the production of a common assessment tool for the European eel in south-western Europe. This question of confidence in complex models as tools used to support fisheries management is crucial, and there is a need to improve transparency regarding model assumptions and uncertainties as well as a need to improve the transfer of model results [12]. The approach used here addressed all these points. To involve stakeholders, it was necessary to present the model assumptions during different workshops, to clarify how these assumptions can affect the model results and to be trans-parent regarding how the stakeholders opinions can be used in the model. The Bayesian approach allowed for the rigorous quantification of the uncertainties, and comparisons of different types of priors allowed us to explore the effect of considering stakeholders opinions. In the article, we focused on the tuning of the model, but the transfer of the results is also a key issue for complex models [12]. In view of this, the most relevant indicators of the model were also discussed in a dedicated workshop in which a questionnaire was used to collect the expectations of the participants about a user friendly web application to display the results. This result in a first web app version which was shown in the previously mentioned Youtube video (https://www.youtube.com/watch?v=6KlhH6FDmBM, accessed on 19 April 2021). We will only be able to determine in the long term whether our participatory approach has reconstructed trust among actors. However, the coordinator of SUDOANG presented the preliminary results of the GEREM and EDA implementation, and the Spanish Ministry found the results and the overall approach interesting. Moreover, a questionnaire sent to participants just after the presentation of the first results showed that most of the participants, at least partially, understood the basis of the model (12 answers out of 16), were satisfied by the zone definition (12 answers out of 14) and were satisfied with their implications in the model construction (10 out of 13; the 3 others were moderately satisfied). Since most of the "negative" answers were still related to minor misunderstandings about the model assumptions, we produced a final YouTube video titled "GEREM in a nutshell", which is visible online (https://t.co/HTPnW9oNz2?amp=1, accessed on 19 April 2021. Sudoang was based on a large diversity of partners (see supplementary material) with different functions: scientists, control forces, fishermen representatives, NGOs, management. With such a large panel, it was difficult to achieve a perfectly balanced panel of types of actors and nationalities at each step of the project. Moreover, in many situations, it was chosen to keep questionnaires anonymous to facilitate participation and avoid selfcensorship. However, this can cause biases in our analysis since we cannot ensure that the respondents were representative of the project. However, a sign-in sheet was used at each meeting so that we were sure that each type of actor has the opportunity to answer and get feedback on the results. Moreover, it confirms that the response rate was quite good compared to the attendance rate.

DISTATIS and the Construction of Stakeholder-Based Priors: Two Relevant Methods Used to Collect Expert Opinion
The Bayesian framework allows for expert knowledge to be used to build informative priors. Here, we follow the approach of [41] to build priors through the analysis of questionnaires sent to stakeholders. The expressed opinions were heterogeneous and thus resulted in very large priors. This is due to the approach used here, in which we only asked the participants to eliminate impossible values, not to provide precise estimates. While this may result in a loss of information, it allowed for more participants to express their opinions and, consequently, to participate in the model construction.
Exploitation rate is not easy to estimate since, theoretically, it requires estimating the abundance of incoming glass eels, a quantity unavailable in most river basins. In this context, stakeholders could only make an estimate of the efficiency of the fishery given its characteristics. Moreover, the estimates which are available in a few river basins are very heterogeneous: more than 95% in the Vilaine River [51], from 15 to 20% in the Loire estuary and 1 to 30% in the Adour estuary [32].
Eliciting expert knowledge requires knowing (1) how the expert opinion will be used, (2) what quantity to elicit, (3) designing the process-here the questionnaire, (4) performing the elicitation and (5) translating the elicited information into a model. Here, the integration of the spatial structure and priors for exploitation rates has a direct use in model construction, so Steps (1), (2) and (5) are quite straightforward. However, the elicitation could have been let differently, for instance, in the future it may be possible to build more informative priors by trying to improve the consensus among participants, for example, by using a DELPHI approach [52,53], which involves several rounds of assessments and discussion among experts; however, this may lead to the loss of some "shy" stakeholders and the risk of overdominance by other stakeholders. Another way of narrowing down the expert opinion range would be to weight the experts, based on a questionnaire for which the answers are known. This weight could then be used when pooling together the different responses [54,55]. However, since our aim was to increase the commitment of experts, we feared that experts receiving too little weight would reject the approach.
The DISTATIS method [34] is a method that is generally used in sensory and consumer science. Here, this method proved to be very relevant when comparing the classification maps made by different stakeholders and was used to achieve a compromise map that allowed us to objectively subdivide the study area into a set of zones. In this exercise, stakeholders opinions showed a good level of consistency with environmental data. For the few participants who had chosen not to remain anonymous, it was possible to see how their maps were located within the "map of experts". In the future, it will be possible to compare the maps drawn by different types of actors and to explore, for example, whether the perceptions of the study region are consistent among French, Portuguese and Spanish stakeholders.
The same zones were also used in the EDA model [27], which was utilised to assess silver eel escapement in the SUDOANG project. However, while recruitment is mostly affected by large-scale factors, eels are subjected to heterogeneous environmental conditions and anthropogenic pressures [16]. Therefore, it would be worthwhile in the future to renew the exercise by asking stakeholders to make new classifications based on not only on their perception of environmental factors in continental waters but also on their perception of human factors: pressures and administrative zoning.
Regarding more specifically the replication potential of the approach with the model GEREM, GEREM was already applied at the French level [28] and on the European scale [21]. In these two analyses, uninformative priors were used for scaling factors and expert knowledge was not accounted for in the zone definition. The possibility to use informative priors on the scaling possible factor makes the model even more flexible as the need for absolute estimates of recruitment in certain river basins becomes less critical: commercial catch data or trap counting are much more frequent than absolute estimates. In this context, the possibility to use the model for the American eel has been discussed [56] and first intents have been made to use GEREM in Japan for the Japanese eel (Yokouchi et al., unpublished).

The Effect of Incorporating Stakeholder Knowledge
Since most of the time series used in this exercise were already used by [21] and are annually used by [23], the resulting trend was very similar to those of previous studies and confirms a decline since the early 1980s. While stakeholders and environmental data suggested to split the zone ATL_IB into two zones (supplementary material Figure S4), it was not possible because of the lack of data in Southwestern Iberian Peninsula. This, combined with the large credibility intervals in some zones (including ATL_F in recent years), urges the implementation of new monitoring series. In this sense, the coordinated monitoring network implemented in the SUDOANG project appears especially relevant and may be used to fill this gap and provide consistent and comparable time series of abundance and biometry. This coordination is especially important since it is well known that inconsistencies in protocols can lead to biases in trend assessment [57].
The results obtained when using stakeholder-based priors were not very different from the results obtained with "control" priors. Obviously, this is due to the wideness of these priors, but this result also shows that stakeholders' opinions were consistent with the data observations. The overall area trend and zone-scale trends were similar regardless of the type of priors used, although recruitment was slightly lower with the "control" priors. At the river basin level, this result also led to differences in the posterior distributions of exploitation rates. While it is unsurprising that the use of stakeholder-based priors did not change the posterior distribution of exploitation rates in river basins for which absolute time series of recruitment are also available (Gironde estuary and Adour River), differences were observed in some other rivers (e.g., the Minho River, Albufera de Valencia lagoon). The two procedures for constructing stakeholder-based priors led to similar priors ( Figure 6) and consequently had very limited impacts on the results. However, we chose to remain very cautious in the construction of the priors. Some authors have recommended restricting the number of experts involved in the construction of priors to find a compromise between redundant information and a good representation of opinion inter-variability [58,59]. By limiting the number of stakeholders, we might have achieved more informative priors and would probably have obtained greater differences between the "triangular" and the "uniform" priors; however, this would require a non-anonymous questionnaire to obtain sufficient information with which to operate stakeholder selection.

Towards the Comparison of Recruitment and Escapement
The SUDOANG project will provide estimates of both recruitment through the use of GEREM [28] and escapement with EDA [27]. By doing so, managers will have comparable and consistent indicators of population statuses over the south-western European area. This is a step forward in fulfilling the need for the coordination of local assessments and to orchestrate initiatives among countries. Currently, we are still not able to carry out the exercise over the whole distribution area [20], but a roadmap has been established to achieve that goal [60].
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/w13091136/s1, Figure S1: Example of map drawn by an expert during the kick-off meeting in Lisbon, Figure S2: Environmental conditions over the study area, Figure S3: Environmental conditions characterizing river basins, Figure S4: The Distatis analysis of maps produced by experts, Table S1: List of partners, An example of questionnaires on exploitation rate is also provided. Funding: The SUDOANG project was funded by the Interreg Program Sudoe V (project SOE2/P5/E0617) through the European Regional Development Fund (ERDF).

Institutional Review Board Statement:
No new monitoring data were collected for the purpose of this study.

Informed Consent Statement:
Informed consent was obtained from all subjects involved in the study.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author. The data are not publicly available due to privacy of some data.