Amandine Courte

Nathalie Cialdella

Alexandre Muller

Vincent Blanfort

Jean-Luc Bochu

Michel Brossard
email: michel.brossard@ird.fr

Identification and evaluation of agricultural practices that store soil carbon, a first step towards low impact agriculture

Keywords: système de culture, Amazonie, petite agriculture tropicale agricultural system, Amazonia, small tropical agriculture

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Le développement de l'agriculture en zone intertropicale concerne 130 000 km 2 .an À1 et engendre des émissions de CO 2 suite aux pratiques de la défriche forestière [START_REF] Giec | Climate Change and Land: an IPCC special report on climate change, desertification, land degradation, sustainable land management, food security, and greenhouse gas fluxes in terrestrial ecosystems[END_REF], tout en étant atteint par les effets du changement climatique (CC). Le sol est un système qui tient une place centrale dans la régulation des écosystèmes agricoles et forestiers (Baveye et al., 2016 ;FAO, 2018), du changement climatique [START_REF] Lal | Soil carbon sequestration to mitigate climate change and advance food security[END_REF] et de la sécurité alimentaire [START_REF] Mcbratney | The dimensions of soil security[END_REF]. En étant à la fois un lieu de stockage du carbone (C) et une source d'émission de gaz à effet de serre (GES), les sols sont l'interface de phénomènes accentués lors des changements d'usage, défriches et mises en culture. Ces modifications accélèrent l'oxydation de la matière organique, et donc altèrent le potentiel de fertilité senso latu des sols (FAO, 2015).

En Guyane, les changements d'usage des terres représentent 43 % des émissions de GES du territoire (GEC, 2016). Les surfaces (< 10 ha par producteur) ne suffisent pas à alimenter la population ; une grande partie des produits alimentaires sont importés ; de plus, l'agriculture est peu rémunératrice [START_REF] Iedom | Rapport annuel 2018 de l'IEDOM[END_REF]. L'agriculture guyanaise doit se développer afin d'assurer l'alimentation d'une population croissante, qui a doublé en 20 ans (INSEE, 2018). Elle doit garantir un meilleur niveau de vie pour les familles, faire face à des contraintes agronomiques et environnementales fortes [START_REF] Renoux | Recenser et évaluer les pratiques agricoles qui stockent le carbone des sols, premier pas vers une agriculture à faible impact en Guyane[END_REF], tout en affichant un impact carbone le plus réduit possible. Une agriculture moins traditionnelle s'est développée dans les années 1960 avec le Plan vert [START_REF] Piantoni | La question migratoire en Guyane française[END_REF], mais reste encore méconnue. Le développement de l'agriculture entraîne sur ce territoire des changements d'usage des terres qui induisent des pertes de C dans la biomasse, mais aussi dans le sol. Dans ce cadre, les agriculteurs doivent mobiliser des pratiques qui puissent maintenir, voire reconstituer, les niveaux du stock organique de carbone (SOC). Ce sont ces pratiques que nous proposons de désigner comme à faible impact sur le carbone et qui nous intéressent dans cette étude.

À l'instar de réflexions similaires sous d'autres latitudes [START_REF] Razafimbelo | Impact de l'agriculture climato-intelligente sur les stocks de carbone organique du sol à Madagascar[END_REF], l'objectif est de s'intéresser aux pratiques à faible impact sur le carbone, dans une réflexion systémique qui prenne en compte la connaissance des systèmes de production. Aussi, nous avons caractérisé les systèmes de production (SP) de la petite agriculture, les pratiques usitées par les agriculteurs et la place du sol dans les prises de décision des agriculteurs et dans la dynamique de gestion des exploitations. Cela permet de savoir si les pratiques à faible impact existent, d'en faire un premier inventaire et de proposer un référentiel. Nous proposons ici une identification et une évaluation de pratiques agricoles vertueuses en termes d'impact sur le climat par le biais d'une vision systémique. À partir de cette analyse des critères agraires, nous avons sélectionné des sites représentatifs pour évaluer les stocks organiques de carbone dans cette petite agriculture. Des mesures de stock ont donc été réalisées de manière synchrone sur ces sites forestiers qui constituent des T 0 (situation à l'état initial) et sur des parcelles agricoles ; ainsi, nous discutons ces pratiques et l'état des stocks de C agricole pour déterminer le devenir du carbone après mise en culture.

Démarche méthodologique 2.1 Les enquêtes

Les enquêtes ont été faites auprès de 69 agriculteurs (soit 1,2 % des exploitants de Guyane) choisis aléatoirement le long du littoral (Fig. 1) et ont permis de décrire les systèmes de production et les pratiques. Les exploitations se situent dans des zones désignées par les services de l'Etat comme « agricoles prioritaires » ou « occupées ». Les premières (prioritaires) sont vouées à être défrichées à l'horizon 2030 pour y développer l'agriculture tandis que les secondes (occupées) sont des zones où se concentre déjà l'activité agricole. Ces exploitations se localisent au Nord de la Guyane, où se concentrent les activités humaines sur 16 sous-types de sols répertoriés dans la carte pédologique. Nous nous focalisons sur l'étude des systèmes maraîchers et arboricoles (séparés ou mélangés au sein d'une exploitation), qui représentent 70 % des systèmes de production guyanais ; les systèmes pâturés sont quant à eux mieux connus et traités dans d'autres programmes.

Nous nous sommes placés dans une dynamique de connaissance de la diversité des systèmes de production pour identifier les modes d'actions efficaces (Cochet et Devienne, 2006) conciliant préservation environnementale et fertilité des sols. En effet, l'objectif est de mettre en avant des pratiques locales et d'en illustrer les impacts sur le stockage du carbone. Les entretiens semi-directifs [START_REF] Imbert | L'entretien semi-directif : à la frontière de la santé publique et de l'anthropologie[END_REF]Olivier de Sardan, 1995a) ont été conduits à partir d'une grille d'entretien construite sur la base de connaissances expertes, afin d'obtenir des données sur : i) l'organisation de la structure agricole, ii) le parcours de vie de l'exploitant, iii) les aspects économiques, iv) les pratiques agricoles. L'arrêt de l'échantillonnage s'est fait une fois arrivé à saturation des données (Olivier de Sardan, 1995b). La typologie a été construite à dire d'experts [START_REF] Landais | Modelling farm diversity: new approaches to typology building in France[END_REF].

La place du sol dans la prise de décision des exploitants

Nous sommes retournés chez 45 des 69 agriculteurs interrogés pour obtenir des compléments de données sur la place que prend le sol dans leurs prises de décision et sur leur niveau de connaissance de l'objet sol.

La perception du sol a été traitée sur la base des questions suivantes : Quels sont les principaux freins à la gestion de votre exploitation ? Qu'est-ce qu'un bon sol selon vous ? Quels indicateurs vous permettent de mesurer ou d'évaluer la qualité de votre sol ? Quelles pratiques vous semblent les plus efficaces pour améliorer ou conserver la fertilité du sol ? Quelles sont les contraintes à leur mise en place ?

Nous avons cherché à déterminer la place que joue le sol dans la logique d'action de l'agriculteur. L'occurrence des réponses a été traitée et a permis de catégoriser les freins au développement agricole, caractériser leur nature et les classer selon trois niveaux d'importance (agriculteurs faiblement impactés, impactés et fortement impactés). L'absence de réponse à une ou plusieurs questions a aussi été prise en compte lorsqu'elle n'exprimait pas une volonté de ne pas participer à l'étude, mais prenait les formes suivantes : « cela ne m'intéresse pas », « je ne me sens pas concerné » ou « je ne sais pas ».

Évaluer le carbone et l'effet des pratiques

Au total, 31 parcelles agricoles (Fig. 1), situées chez 10 agriculteurs représentatifs de la typologie (en termes de pratique et de gestion) et dans une station expérimentale, ont été sélectionnées pour réaliser des prélèvements de sol (5 parcelles avec usage de pratiques à fort intrants, 20 avec usage de pratiques à faible impact et 6 avec usage de pratiques traditionnelles). Plusieurs parcelles ont été échantillonnées par agriculteur lorsque les sols ou les pratiques variaient au sein d'une même exploitation. De nombreuses parcelles ont été étudiées dans des systèmes à faible impact, puisque ces pratiques nous intéressent particulièrement dans cette étude et que ce sont les pratiques les plus diversifiées. Les mesures du stock organique de carbone agricole et forestier sont classiques (FAO, 2019 ; RMQS : www.gissol.fr) et s'appuient sur des mesures de densité apparente. Ici, les mesures sont réalisées sur des échantillons composites pour les couches 0-0,3 et 0,3-0,5 m de sol, en prenant en compte, le cas échéant, l'organisation spatiale des parcelles (rang et inter-rang). Les comparaisons entre les SOC à T 0 avant défriche et les SOC agricoles à T þn (entre 2 et 30 ans) ont été effectuées à masse de sol équivalente pour la couche considérée. Les échantillons T 0 ont été choisis à dire d'experts en fonction des particularités pédologiques observées sur le terrain et prélevés au voisinage des sites. Nous étudions donc des pratiques mises en place à des pas de temps différents. Le bilan des émissions des gaz à effet de serre a été effectué chez 9 de ces 10 agriculteurs en s'appuyant sur l'outil ACT-DOM. L'outil est en cours d'adaptation sur le territoire guyanais et permet d'établir des diagnostics Energie et GES à l'échelle de la ferme (Dallaporta et al., 2016). Pour utiliser cet outil, nous avons collecté des données quantitatives (intrants, surfaces couvertes, densité du plastique des serres, quantités produites et exportées, etc.). Ainsi, nous obtenons pour chaque exploitation ses émissions en teqCO 2 .ha À1 .an À1 . Ce bilan n'a pas pu être réalisé chez des agriculteurs aux pratiques traditionnelles (abattis-brûlis itinérants), car les données quantitatives sont absentes. Cela constitue une limite pour l'obtention de bilans carbone fins.

Résultats

Les systèmes de production observés

Les systèmes de production observés se distinguent dans un premier temps par les pratiques (à faible impact, à forts intrants et traditionnelles) et dans un second par l'objectif de production (commercialisation, transformation, vivrier) (Fig. 2).

Les systèmes à faible impact (30 exploitations)

Dans ce groupe, nous retrouvons des exploitations récentes (maraîchage et arboriculture) en phase d'installation (24 %) et des exploitations plus anciennes issues d'un héritage, où le système de production est bien en place. Le projet partagé de ces exploitants est d'optimiser leur production dans une démarche de respect de l'environnement ; 30 % d'entre eux sont d'anciens agriculteurs conventionnels (faisant usage de pratiques à forts intrants) qui ont réalisé une conversion vers l'agriculture biologique. La prise de décision est fonction de valeurs et n'est pas dictée directement par un objectif d'augmentation du chiffre d'affaire ou de la production. Les agriculteurs y développent de nouvelles pratiques avec une gestion plus organique et sans produits chimiques. Ils travaillent essentiellement manuellement et réutilisent les sous-produits de leur exploitation ; 10 % produisent et utilisent du bois raméal fragmenté (BRF) issu de la défriche de leur sous-bois. Si les ressources financières pour engager un salarié le permettent, l'entretien des parcelles est alors partiellement délégué. Lors des périodes de récolte, la famille ou les amis viennent soutenir l'activité. Leur production de fruits ou légumes est vendue sur les marchés, directement depuis l'exploitation, à la coopérative locale ou dans de grandes et moyennes surfaces. Sont comprises dans ces systèmes les exploitations en polyculture-élevage et des transformateurs (commercialisation d'un unique produit transformé) ; l'activité d'élevage ou de transformation y représente une grande partie du revenu, du temps ou de l'investissement. En polycultureélevage, l'atelier animal soutient économiquement l'atelier végétal. Les parcelles sont organisées dans l'espace et les pratiques orientées vers la fertilité du sol : utilisation de fumier, couvert végétal, rotation, bois raméal fragmenté. Dans les systèmes de production dédiés à la transformation, la production des produits transformés représente l'intégralité du revenu. Lorsque la production ne suffit pas à alimenter l'atelier de transformation, l'approvisionnement est complémenté par des achats extérieurs. Dans les systèmes de production avec un objectif d'autoconsommation, les agriculteurs bénéficient souvent d'un revenu supplémentaire à la simple vente des surplus de leur production, leur permettant d'investir dans leur activité agricole. Les pratiques agricoles sont réfléchies, se veulent à faible impact et particulièrement diverses et les caractéristiques du sol sont maîtrisées. Ces agriculteurs n'ont pas de grandes surfaces (moyenne de 8,5 ha) et leurs parcelles se situent à proximité de l'habitation.

Les systèmes à forts intrants (19 exploitations)

Dans ces exploitations mécanisées, le recours aux intrants chimiques est systématique et s'accompagne rarement d'apport de matière organique fraîche. Le temps de travail est partagé entre maraîchage et arboriculture, les parcelles sont généralement monospécifiques (ananas, concombre, banane, agrumes, etc.). Le revenu est exclusivement issu de l'activité de l'exploitation et est suffisant pour lui permettre de se développer.

Les systèmes traditionnels (20 exploitations)

La particularité de ces systèmes de production réside dans la difficulté d'accès au foncier. Le premier groupe d'exploitations (40 %) correspond à des zones d'habitation traditionnelle très ancienne, mais sans titre de propriété ni statut agricole. L'agriculture est plutôt vivrière et ne représente pas l'activité principale. Ces exploitants n'ont pas les moyens d'investir dans l'agriculture. Le second groupe (60 %) qui se distingue dans ces systèmes de production est formé par des agriculteurs nouvellement installés qui perçoivent les aides et subventions sectorielles et qui cherchent à diversifier leurs productions, dans des circuits courts ou sur les marchés locaux et grâce à un statut professionnel agricole reconnu. Les exploitations se situent sur des terrains contraignants, aux sols sableux ou aux pentes fortes. Les surfaces cultivées sont petites (< 2 ha) et la production manuelle plutôt diversifiée : tubercules et fruits. Dans la moitié des cas, il s'agit d'abattis-brûlis avec une tendance à la diminution des surfaces cultivées. L'utilisation d'intrants est occasionnelle en fonction de la trésorerie disponible. Des planches de maraîchage, constituées d'un mélange de terre et de compost, leur permettent de pallier les contraintes des sols peu fertiles. Les techniques utilisées proviennent des acquis familiaux.

Les freins rencontrés par les systèmes de production

Nos résultats montrent que les agriculteurs les plus efficients sur le plan économique, aussi bien dans les « systèmes de production à forts intrants » que dans les « systèmes de production à faible impact », arrivent à organiser, et parfois à gérer eux-mêmes, la commercialisation de leur production grâce à une main-d'oeuvre familiale importante dans le cas des « SP à forts intrants », ou grâce au recours à des marchés de niche pour les « SP à faible impact ». Plusieurs freins à l'agriculture et au développement de pratiques à faible impact sont mis en évidence (Tab. 1). La lourdeur des démarches administratives, le manque d'instituts techniques et de filières structurées est manifeste. L'enclavement géographique, les terrains isolés difficiles d'accès compliquent la commercialisation ; or la vente sur les marchés urbains constitue une des seules options rentables. Les conditions pédoclimatiques contraignantes sont systématiquement citées parmi les difficultés, même si elles sont moins importantes dans les systèmes de production à forts intrants. La saison des pluies, avec des précipitations de 2000 à 4000 mm. an À1 , est responsable d'excès d'eau et de maladies, et l'intensité des pluies favorise le ruissellement, puis l'érosion des sols. La couche supérieure du sol où est localisée la matière organique est en général peu épaisse, les sols sont naturellement acides, souvent carencés en calcium, ce qui nécessite des amendements calciques réguliers et coûteux. La variabilité des sols sur une même parcelle rend également difficile la gestion des cultures.

Description de la diversité de pratiques

On distingue un large panel de pratiques (Tab. 2). Les plus répandues, présentes chez 42 % des agriculteurs interrogés, que nous nommerons « à forts intrants » par commodité, consistent en l'usage d'intrants chimiques pour la gestion des ravageurs et des maladies, des adventices et pour la fertilisation. Ces pratiques sont mobilisées dans les systèmes de production à forts intrants et ne sont pas détaillées ici, celles qui nous intéressent le plus étant celles à faible impact.

Sur la population étudiée, 51 % des agriculteurs mettent en place au moins une pratique de conservation ayant pour but d'augmenter la fertilité du sol et limiter les émissions de GES. Ces pratiques, décrites ci-après, peuvent-être associées pour un effet combiné : utilisation des débris de fauche et de paillis ( 20

À forts intrants þ þ þþþ þ þþ À faible impact þþ þ þþ þþ þþ Traditionnel, vivrier þþþ þ - þ - À Pas impacté. þ faiblement impacté. þþ Impacté. þþþ Fortement impacté.
jachère de durée variable (27 % des agriculteurs). La parcelle est mise en repos pendant 2 ou 3 ans, voire pendant plus de 7 ans, pour permettre la croissance du recru forestier, qui est ensuite brûlé avant de réinvestir la parcelle ; confection d'un compost amélioré (5 % des agriculteurs).

Un mélange de charbon, débris végétaux et fumier est constitué pour faire une « terre noire » inspirée des sols amérindiens précolombiens (terra preta). Certains peuvent aller jusqu'à la réalisation de charbon artisanal pour confectionner ce mélange, ce qui questionne le bilan GES d'une telle pratique.

La perception des sols par les agriculteurs

Les agriculteurs observent principalement leurs sols à travers la texture. Ils distinguent : « les sols rouges argileux », « les sols argileux mais un peu plus sableux », « les sols à dominante de sables » (ferralitisols argileux et argilo-sableux) et « les sols à sables blancs » (fluvisols, podzosols). Le drainage est noté dans la différenciation des sols, ce qui va déterminer l'emplacement des différentes cultures : « sols argileux drainants et trop secs en saison sèche », « sols hydromorphes » ou « gorgés d'eau en saison des pluies ». La présence d'éléments grossiers est souvent citée comme élément favorisant le drainage, même si aucune différence n'est évoquée quant à la nature et la taille de ces éléments (cailloux et gravillons ferrugineux, quartz...). Le travail du sol étant souvent concentré sur de petites surfaces peu mécanisées, graviers et cailloux sont souvent considérés comme un atout plutôt qu'un frein à la mise en culture. La couleur du sol est systématiquement observée, avec l'idée que « plus la terre est sombre, plus elle est riche en matière organique ». La structure grumeleuse est citée comme un paramètre indiquant la bonne qualité des sols, car cela démontre une bonne activité biologique et un sol « qui respire » (ferralitisols à drainage vertical libre).

Les agriculteurs ayant suivi des formations prêtent souvent plus attention aux indicateurs chimiques des sols, les plus importants et les mieux maîtrisés étant la matière organique (MO) et le pH. L'amendement calcique est réalisé par 92 % des agriculteurs. Seuls 7 % ne gèrent pas la fertilité de manière chimique ou organique (Tab. 2).

Enfin, la profondeur utile du sol n'est pas considérée comme un élément important. Les agriculteurs estiment que leurs sols sont profonds (sauf présence de nappe dans les sols hydromorphes) et homogènes.

Ces perceptions générales et l'intérêt porté à ces variables par les agriculteurs varient pour chaque groupe, en fonction notamment de leur objectif de production.

Dans le cadre des systèmes de production à forts intrants, les agriculteurs considèrent qu'un « bon » sol est productif, leur indicateur est donc le rendement qui est un critère subjectif la plupart du temps car il n'est pas mesuré. La pratique avancée comme solution à leurs sols peu fertiles est l'utilisation de fertilisant. Ils se déclarent désintéressés par la question de la qualité des sols. Le sol est alors plus un support qu'une ressource. C'est d'ailleurs ce groupe qui considère les conditions pédoclimatiques comme un frein léger à la mise en culture (Tab. 1).

Dans le cadre des systèmes de production à faible impact, les agriculteurs estiment qu'un « bon » sol est productif et le voient comme un objet vivant nécessaire à comprendre. Leurs indicateurs sont variés et souvent liés à l'aspect général du sol (texture, couleur, pierrosité, parfois structure, conditions hydrodynamiques [présence ou absence de phases inondées ou à l'inverse de sécheresses, capacités drainantes]). Certains citent le pH et les analyses de sol. Ils font appel à des pratiques variées et combinées.

Dans le cadre des systèmes de production traditionnels, les agriculteurs observent la quantité et la qualité des produits récoltés pour déterminer si le sol est « bon ». Ils se basent tous sur la couleur et la texture pour décrire leur sol et parfois sur la quantité d'adventices qui déterminent le moment où il va falloir changer de parcelle dans une dynamique d'abattis itinérant. Ils disent mobiliser un ensemble de connaissances transmises par leurs ancêtres (buttes, paillis, brûlis, jachère) pour garantir la qualité de leurs terres sur le long terme.

Tableau 2. Occurrence (%) des réponses des exploitants pour les différents modes de gestion et pratiques agricoles. Le labour est limité aux 30 premiers centimètres de profondeur. Dans les enquêtes, le facteur sol est apparu déterminant. Dans la gestion spatiale du système de production, nous observons que l'agriculteur s'adapte aux spécificités des sols par le choix des cultures et des pratiques. Par ce biais, l'objectif est bien entendu d'optimiser la production en bénéficiant des « potentialités naturelles » du sol. Ainsi, concernant la couleur du sol, il apparaît que 90 % des agriculteurs y sont sensibles et mettent donc indirectement la teneur en MO en avant, clé du maintien de la fertilité des sols.

Effet des pratiques sur le carbone du sol

Pour la couche 0-0,5 m, on observe une grande diversité de valeurs du stock organique de carbone dans les parcelles agricoles (Tab. 3), de 30 à 425 t C.ha À1 , tout comme leurs références sous forêt, où les SOC sont compris entre 51 et 371 t C.ha À1 . Les SOC agricoles les plus élevés correspondent à des sols cultivés sur des sols aux SOC naturellement élevés sous forêt, dont la défriche récente est associée à l'utilisation de bois raméal fragmenté (paillis constitué grâce au broyage des bois de défriche). A contrario, les déstockages les plus importants sont représentés par des parcelles en maraîchage avec utilisation d'intrants chimiques ou par l'arboriculture à faible impact. Cela est lié à deux facteurs : la nature sableuse des sols (podzosols et fluvisols) et des teneurs en matière organique et éléments minéraux assez réduites ou un apport de matière organique limité. En effet, la diminution des entrées de carbone du sol est ici également déterminante, comme observé par [START_REF] Fujisaki | Soil carbon stock changes in tropical croplands are mainly driven by carbon inputs: A synthesis[END_REF] sur la station expérimentale de Combi.

Dans les systèmes de production à faible impact, la tendance est à l'accroissement des stocks avec le temps. En effet, lorsque la défriche et les pratiques sont à faible impact, les parcelles étudiées se rapprochent de leur niveau de stock initial, même à très court terme. Sur les 10 systèmes maraîchers à faible impact étudiés, les 6 qui utilisent le bois raméal fragmenté ont dépassé le niveau de stock initial. On peut alors se poser la question des transferts de fertilité par les amendements organiques réalisés. Ces résultats doivent être mis en perspective avec les autres travaux amazoniens aux variations de stocks hétérogènes allant de þ20 à À40 % après mise en culture (Fujisaki et al., 2015a). Il existe encore trop peu de données pour en déduire une règle stricte en dehors des systèmes pérennes à pâturages pour le bassin amazonien.

Les premiers bilans de gaz à effet de serre des exploitations montrent, en fonction des variations du stock organique de carbone, que les systèmes de production à forts intrants ont plutôt une tendance au déstockage de carbone et sont la forme d'agriculture la plus émettrice (Tab. 3). Les résultats pour l'arboriculture sont variables, mais les bilans d'émission sont très faibles. S'agissant du maraîchage à faible impact, les émissions se situent entre 1 et 3 teq CO 2 .ha À 1.an À1 avec une grande variabilité du stock organique de carbone (liée au nombre d'années de culture et à l'effet de la valeur du stock de carbone de référence initial). Dans l'ensemble, l'usage des buttes ou de terrasses permet un accroissement considérable des SOC entre 12 et 62 % (0-0,5 m de sol), et ce en un laps de temps très court (moins de 5 ans). Cette observation est encourageante, car elle entre dans le principe du 4 ‰ (Derrien et al., 2016), rappelant que l'agriculture n'est pas forcément émettrice de gaz à effet de serre, mais peut au contraire participer à la séquestration du carbone. Toutefois, ces effets sur le stock ne sont pas pérennes du fait de la disparition des fragments de bois, ce qui est inéluctable dans ces conditions climatiques [START_REF] Perrin | Conversion of forest to agriculture in Amazonia with the chop-and-mulch method: Does it improve the soil carbon stock? Agriculture[END_REF]Fujisaki et al., 2015b).

Conclusion

Nous avons identifié 3 systèmes de production caractérisés par des objectifs de production et de pratiques qui mettent en valeur quelques originalités. L'analyse des freins au développement de l'activité agricole et celle de la perception des sols par les agriculteurs montrent que « l'objet sol » est un indicateur reconnu qui conditionne les modes de gestion de l'espace.

Les stocks de carbone organique des sols sous végétation naturelle sont élevés et varient rapidement et fortement après la défriche et la mise en culture. Les mesures chez les agriculteurs sont cohérentes et soulignent bien la fragilité des systèmes solplante dès lors que l'on défriche. Nous avons vu en particulier que les phénomènes sont rapides, tant sur la réactivité à constituer des stocks de carbone exemplaires avec certaines pratiques, qu'à les voir se dégrader dans certains systèmes intensifs. Si l'impact des pratiques à forts intrants se traduit par des baisses drastiques du stock organique de carbone, jusqu'à À 50 % du SOC à T 0 en 5 ans, les pratiques à faible impact peuvent avoir des effets bénéfiques sur les dynamiques de stockage du carbone, en plus d'avoir un impact environnemental réduit.

Nous soulignons en particulier le caractère exploratoire des bilans des émissions des exploitations que nous avons effectués. Ces premiers résultats montrent que les émissions sont les plus élevées dans les sites à forts intrants. Notre étude montre qu'il est possible d'avancer dans le secteur des terres entre carbone environnemental et carbone agricole et qu'il est possible de mettre en place des solutions durables. La Guyane s'inscrit aussi bien dans les préoccupations des pays en développement que du Nord, qui insistent depuis la Cop 17 de 2011 sur l'adaptation et l'atténuation. Quels que soient les choix de transformation de la forêt et des sols, cela reste un enjeu.

Références

Baveye PC, Baveye J, Gowdy J.

Fig. 1 . 1 .

 11 Fig. 1. Situation géographique des exploitations étudiées et des sites où les stocks de carbone des sols ont été prélevés en milieu agricole. Fig. 1. Map of farms surveyed and agricultural places where soil carbon stocks have been measured.

Fig. 2 .

 2 Fig. 2. Typologie des exploitations. Fig. 2. Farm typology.

Table 1 .

 1 Tableau 1. Hiérarchisation des freins au développement agricole. Impediments to farming development.

	% des

tration de la matière organique. Ces buttes sont généralement décalées dans l'inter-rang tous les trois ans pour garder une parcelle relativement fertile ;

Table 2 .

 2 Occurrence (%) of farmers' responses for different modes of management and farming practices.Tableau 3. Stocks de carbone (t C.ha À1 ) des sols agricoles et variations par rapport à la référence sous forêt (n = 31 parcelles), couche 0-0,5 m, en fonction du temps et du mode de culture.

	Gestion des ravageurs/maladies	Aucune	43
		Organique	34
		Chimique	23
		Aucune	7
	Gestion de la fertilité	Organique	51
		Chimique	42
		Aucune	13
	Gestion des adventices	Organique	57
		Chimique	30
		Aucune	73
	Jachère	0-7 ans	18
		þ de 7 ans	9
	Irrigation	Oui Non	57 43
	Paillis	Oui Non	20 80
	Plantes de couverture	Oui Non	7 93
		Brûlis þ bulldozer	18
	Mode de défriche *	Brûlis Bulldozer	60 12
		BRF	10
	Utilisation de « Terra preta »	Oui Non	5 95
	Association d'espèces	Oui Non	22 78
	Buttes	Oui Non	9 91
	*Brûlis þ bulldozer : utilisation combinée des deux modes de
	défriche. Le feu permet un éclaircissement du sous-bois avant le
	passage du bulldozer pour coucher les bois restants et procéder au
	débardage.
	Brûlis : Bois entièrement brulés sur la parcelle.
	Bulldozer : Bois couchés et débardés de manière mécanique.
	BRF : le sous-bois et les houppiers sont broyés et constituent un mulch
	avant le passage d'engins légers pour l'exploitation des troncs
	> 20 cm de diamètre.

Table 3 .

 3 Soil carbon stocks (t C.ha À1 ) in farm plots, 0-0.5 m layer, in function of time and practices.

	Stock	Stock	Delta	Delta	GES (teq	Temps depuis mise	Information	Type de sol
	agricole	forestier	Stock	Stock (%)	CO 2 .ha À1 .an À1 )	en culture (année)
	Systèmes à forts intrants
	100,05	102,55	À2,5	À2,44	4,5	20	Reprise de la parcelle après une	Ferallitisol
							jachère de 2 ans pour du maraîchage
							monospécifique, rotations
	36,92	70,34	À33,42 À47,51	6,3	8	Maraîchage sous serre,	Redoxisol
							monospécifique
	37,92	70,34	À32,42 À46,09	6,3	8	Maraîchage sous serre,	Redoxisol
							monospécifique
	29,68	73,02	À43,34 À59,35	nd	30	Maraîchage monospécifique	Thalassosol
	143,21	161,34	À18,13	11,24	nd	5	Arboriculture monospécifique	Ferallitisol
	Systèmes à faible impact
	Maraîchage
	206,85	128,97	77,88	60,39	2,5	5	Culture sous serre, sur terrasse, BRF,	Ferallitisol
							rotations
	123,45	128,97	À5,52	À4,28	2,5	2	Culture sur terrasse, rotations	Ferallitisol
	114,60	121,91	À7,31	À6,00	3	15	Buttes, rotations et associations	Ferallitisol
	204,79	126,45	78,34	61,95	1	4	Buttes, BRF, rotations et associations	Ferallitisol
	90,88	80,82	10,06	12,45	1	3	Buttes, BRF, rotations et associations	Arenosol
	241,4	201,09	40,31	20,05	1	3	Buttes, BRF, rotations et associations	Redoxisol
	101,3	126,45	À25,15 À19,89	1,4	10	Buttes avec paillis récent, rotations	Ferallitisol
	126,18	126,45	À0,27	À0,21	1,4	10	Buttes avec paillis décomposé,	Ferallitisol
							rotations
	414,47	371,31	43,16	11,62	1,7	5	Planches avec des apports conséquents	Redoxisol
							de matière (BRF et fumier), rotations
	425,33	371,31	54,02	14,55	1,7	5	Planches avec des apports conséquents	Redoxisol
							de matière (BRF et compost),
							rotations
	Arboriculture
	125,77	126,45	À0,68	À0,54	3	4	Agro-forêt	Ferallitisol
	82,77	111,26	À28,49 À25,61	2	20	Monospécifique	Ferallitisol
	122,4	121,91	0,49	0,40	1	10	Monospécifique	Ferallitisol
	124,8	121,91	2,89	2,37	3	20	Monospécifique	Ferallitisol
	69,82	112,00	À42,8	À37,66	0,3	6	Inter-rang, Monospécifique	Ferallitisol
	123,79	112,00	11,79	10,53	0,3	6	Rang, Monospécifique	Ferallitisol
	56,9	71,61	À14,71 À20,54	0,3	2	Monospécifique	Podzosol
	82,88	51,32	31,56	61,50	nd	30	Agro-forêt	Redoxisol
	Station expérimentale : Combi, grandes cultures et prairies (Fujisaki et al., 2014)
	77,7	73,2	4,5	6,15	nd	5	Prairie	Ferallitisol
	70,5	73,2	À2,7	À3,69	nd	5	Culture de céréales	Ferallitisol
	64,1	73,2	À9,1	À12,43	nd	5	Semi direct de céréales	Ferallitisol
	Systèmes traditionnels
	93,47	94,00	À0,53	À0,56	nd	2	Changement de parcelle tous les 3 ans	Ferallitisol
	51,95	78,91	À26,96 À34,17	nd	30	Abattis	Ferallitisol
	61,84	65,26	À3,42	À5,24	nd	30	Abattis	Ferallitisol
	93,69	87,07	6,62	7,60	nd	30	Abattis	Ferallitisol
	63,56	73,21	À9,65 À13,18	nd	30	Abattis	Ferallitisol
	68,31	60,25	8,06	13,38	nd	15	Abattis	Ferallitisol

  2016. Soil "Ecosystem" Services and Natural Capital: Critical Appraisal of Research on Uncertain Ground. Front. Environ. Sci. 4. DOI: 10.3389/fenvs.2016.00041. Cochet H, Devienne S. 2006. Fonctionnement et performances économiques des systèmes de production agricole : une démarche à l'échelle régionale. Cahiers Agricultures 15: 578-583. DOI: 10.1684/agr.2006.0028. Dallaporta D, Bochu JL, Vigne M, Ouliac B, Zoogones LJ, Lecomte P, et al. 2016. Taking into account carbon sequestration of pasture in carbon balance of cattle ranching systems established after deforestation in Amazonia. In: Proc. 10th Int. Rangeland Congress, Saskatoon, Canada, 6-22 July 2016, pp. 399-400. Derrien D, Dignac MF, Basile-Doelsch I, Barot S, Cécillon L, Chenu C, et al. 2016. Stocker du C dans les sols : quels mécanismes,

A.Courte et al. : Cah. Agric. 2020, 29, 21


Gilles Pagès

Olivier Pironneau

Stochastic Control: from Gradient Methods and Dynamic Programming to Statistical Learning

Keywords: Stochastic control, Dynamic Programming, Partial Differential Equations, Statistical Learning

In this article the authors wish to contribute to the validation of supervised learning for stochastic control. We review the popular methods for stochastic control and compare their numerical performance to those of a neural network. This will be done on a simple but practical example: optimal fishing quotas to preserve the biomass of a fishing site.

Introduction

Stochastic control of differential systems with random data or modelling errors, is usually solved by Dynamic Programming, a technique introduced in the 1950 by Richard Bellman [START_REF] Bellman | Dynamic Programming[END_REF]. Numerically one has to solve a backward in time differential equation possibly coupled with the forward stochastic differential equation (SDE) of the problem. The method has various discrete implementations (see [START_REF] Pagès | An Optimal Markovian Quantization Algorithm For Multi-Dimensional Stochastic Control Problems[END_REF], [START_REF] Lauriere | Dynamic Programming for mean-field type control[END_REF], etc). Calculus of variations and optimal control is another approach, which, in principle is more expensive because the backward adjoint Partial Differential Equation (PDE) and the forward Kolmogorov PDE for the PDF of the process are coupled and must be solved iteratively by a gradient method. A third method has been proposed several years ago (see [START_REF] Balakrishnan | Adaptive-critic-based neural networks for aircraft optimal control[END_REF], [START_REF] Kamalapurkar | Reinforcement Learning for Optimal Feedback Control[END_REF], [START_REF] Bertsekas | Reinforced Learning and Optimal Control[END_REF], [START_REF] Auger | Parameter Identification by Statistical Learning of a Stochastic Dynamical System Modelling a Fishery with Price Variation[END_REF] and others) : Supervised Learning with Neural Networks: like off-line/on-line reduction methods, one is ready to pay a high CPU price in a "training" phase so as to solve quickly the problem with new data with the "trained network". The novelty is that a trained neural network can be ported on a small unit like a smart phone and yet be fast and accurate. In this article the three classes of methods will be tested on an easy yet practical problem studied in [START_REF] Brochier | Can overexploited fisheries recover by self-organization? Reallocation of the fishing effort as an emergent form of governance[END_REF]: the control of a fishing site in Senegal where the fish biomass is threatened by intensive fishing. In [START_REF] Auger | Parameter Identification by Statistical Learning of a Stochastic Dynamical System Modelling a Fishery with Price Variation[END_REF] it was shown that parameter identification from a few measurements near time zero could be obtained accurately and fast by supervised learning. Here the parameters are known but with some uncertainty and the objective is to apply optimal quotas to stabilize the fish biomass at a desired level in this stochastic setting. The model has two ordinary differential equations (ODE), one for the fish biomass X t at time t, and one for the fishing effort which, for us, will be the number of boats at sea E t . 1 gilles.pages@sorbonne-universite.fr, LPMA, Sorbonne Université, Paris, France 2 olivier.pironneau@sorbonne-universite.fr , LJLL, Sorbonne Université, Paris, France.

Quotas are meant to keep the biomass at a constant desired level X d . Hence the total number of fish caught must be lower than the quota Q t . By construction each boat can't catch more than qX t , where q is a constant reflecting the "structural limit" of the fisherman's catch yet proportional to the biomass. Nevertheless each boat is told to keep their catch below Q t /E t . The control problem is solved by minimizing a criteria, the integral over time of X t -X d 2 . But quotas, imposed by the regulator, are unpopular. To avoid too strict quotas, we add to the criteria the integral over time of -αQ t where α is constant and positive; consequently, a small quota will not be optimal because it gives a larger value to the criteria than a big quota. Along the way we found that the solution of the problem is very often "bang-bang", i.e; switching from the structural limit to the quota often. For the regulator, rapidly changing quotas may be also dangerously unpopular; hence in an attempt to avoid such frequent drastic changes we add to the criteria the integral over time of β Qt 2 where β is a positive constant. More details on the modelling of the problem can be found in Appendix A.

The plan of the paper is as follows: In section 1, the control problem is stated. It is a stochastic control problem because the data are not known accurately. Existence of solutions and discretization methods are discussed. We identify 3 formulations:

1. Deterministic control where the state is the solution of an ODE. 2. Non-dynamic control, i.e. an SDE for the state equation but a deterministic control function of time only. 3. A fully stochastic state and control, usual to Dynamic Programming.

In section 2, a Monte-Carlo method is proposed for the non-dynamic case. It is a brute force Monte-Carlo approach, much more expensive than the Dynamic Programming approach of section 4, but it has an advantage: it can be ported to a Neural Network for training and testing. We have used the standard multilayered Recurrent Neural Networks of Keras [START_REF] Chollet | Deep learning with Python[END_REF] (RNN, see figure 3); LSTM is also available in Keras, but as there are no repeating patterns in the solution of a SDE with brownian noise, it is not clear that a Long-Short-Term-Memory feedback architecture would improve the (already good) performance of RNN. Finally, in the last section PDE methods are investigated, based either on Ito Calculus or Kolmogorov equations for the probability density function (PDF) of the process.

A model problem

To study the evolution of the biomass of a fishing site we consider the following stochastic optimal control problem

min u∈U J(u) := T 0 E (X t -X d ) 2 -αu(X t , t) + β|∂ t u(X t , t)| 2 dt, (1)
with dX t = X t (r -κX t -u t )dt + X t σdW t , X 0 = X 0 ≥ 0, (2)
where W t is a N 0,1 Weiner process; X d , r, κ, σ, u m , u M are real constants which depend on the fishing site considered, and

U = {u ∈ H 1 (]0, T [) : u m ≤ u(X t , t) ≤ u M , a.e.} (3)
This problem arises for the determination of fishing quotas (see Appendix A for details): X t is the fish biomass at time t and u t X t is the fishing quota at time t. Occasionally we will also display a variable proportional to the number of boats at sea, the "fishing effort" E(t), modelled by

dE t = (a -u t X t -cE t )dt, E 0 = E 0 given.
where a is related to the price of fish and c is the operating cost of the fisherman. As explained in Appendix A, a -u t X t is the profit of the fisherman and it decreases when too many fish are caught because the price of fish collapse.

The control problem aims at keeping a biomass at a desired level X d at all times, but at the same time, quotas being unpopular, there is a term to prevent small quotas, -αu t and another to prevent quotas from changing too quickly: β|∂ t u| 2 . Note that (X, t) → u(X, t) is a deterministic function and {u(X t , t)} t is a Markov process. Finally,

∂ t u(X t , t) means ∂ t u(X, t)| X=Xt .
It is important to notice that one may also consider the case where u is a deterministic function of t only:

min u∈U J(u) := T 0 E (X t -X d ) 2 ] -αu(t) + β|∂ t u(t)| 2 dt, with (2). (4)
We shall refer to it as the non-dynamic problem by opposition to Dynamic Programming where u(X t , t) is stochastic because of X t but Markovian. In the non-dynamic case, the optimal control is the same for all trajectories X t and in the dynamic case u adapts to the trajectories. To our surprise the non-dynamic case turns out to be more precise and more suitable because the day-to-day noise may not be very meaningful.

Occasionally we shall also refer to the deterministic case, i.e. when σ = 0.

Existence of solution: the deterministic case

In the deterministic case the problem is

min u∈U {J(u) = T 0 |X -X d ] 2 -αu + β| u| 2 : Ẋ = X[(r -κX -u), X(0) = X 0 } (5)
As u is positive 0 ≤ X(t) ≤ X M := X 0 e rT . Then obviously F (X) := rX -κX 2 -Xu is uniformly Lipschitz in X:

|F (X 1 (t)) -F (X 2 (t))| ≤ |X 1 (t) -X 2 (t)|(r -u m + 2κX M |, ∀t ∈ [0, T ];
therefore the solution of the ODE exists and is unique. It is also uniformly Lipschitz continuous in u:

Ẋi = rX i -κX 2 i -u i X i , i = 1, 2, ⇒ Ẋ1 -Ẋ2 = (X 1 -X 2 ) (r -κ(X 1 + X 2 ) -(u 1 + u 2 )) -(u 1 -u 2 )(X 1 + X 2 ), ⇒ | Ẋ1 -Ẋ2 | ≤ |X 1 -X 2 |(r + 2X M + 2u M ) + 2X M |u 1 -u 2 |, ⇒ |X 1 -X 2 | ≤ 2|u 1 -u 2 |X M (exp((r + 2X M + 2u M )t) -1)/(r + 2X M + 2u M ) (6)
where the last line is due to the Grönwall lemma. This shows also that u Proof: Consider a minimizing sequence: u n ∈ U, J(u n ) → inf u∈U J(u). U being nonempty, bounded and closed and u 2 0 being part of the criteria, there is a subsequence with u n → u * ∈ U weakly in H 1 (]0, T [) and strongly in L 2 (]0, T [). By the above established continuity, X n → X * in L 2 (]0, T [). Obviously T 0 αu n → T 0 αu * and by the lower semicontinuity of u → u 2 0 ,

→ T 0 |X(t) - X d | 2 dt is continuous in L 2 (]0, T [).
J(u * ) ≤ lim inf J(u n ) = inf u∈U J(u).
Therefore u * is a minimum and problem (5) has a solution. X s ds = +∞ is impossible unless τ = +∞, P a.s. Therefore X t ∈ R + , ∀t and

X 0 e St(u M ) exp -κ t 0 e Ss(u M ) ds ≤ X t ≤ X 0 e St(um) . (8)
Theorem 1. Problem (1)(2) has a solution.

Proof: Consider two trajectories corresponding to the same event but different u, denote δ the difference and ¯the average; we have

dδX t = [δX t (r -2κ Xt -ū) -Xt δu]dt + σδX t dW t .
By Itô calculus:

(δX t ) 2 = 2 t 0 [δX s ) 2 (r + σ 2 2 -2κ Xt -ū) -Xs δu(s)δX s ]ds + 2 t 0 σ(δX s ) 2 dW s . Hence E[(δX t ) 2 ] ≤ 2(r + σ 2 2 -u m ) t 0 E[(δX s ) 2 ] - t 0 δu(s)E[ Xs δX s ]ds ≤ 2(r + σ 2 2 -u m ) t 0 E[(δX s ) 2 ] + t 0 δu(s)E[( Xs ) 2 ]ds. (9)
We have seen earlier that 0 ≤ X t ≤ X 0 e (r-1 2 σ 2 )t+σWt , therefore E[( Xs ) 2 ] ≤ 4e (2r+σ 2 )s E[X 2 0 ], so much so that by Grönwall's lemma,

E[(δX t ) 2 ] ≤ 2(r + σ 2 2 -u m ) t 0 E[(δX s ) 2 ]ds + 4e (2r+σ 2 )T E[X 2 0 ] t 0 (δu(s)) 2 ds;
and by Grönwall's lemma again,

E[(δX t ) 2 ] ≤ 4e 2r+σ 2 -2um)t e (2r+σ 2 )T E[X 2 0 ] t 0 |δu(s)|ds.
Consequently, for some

C ∈ R + function of r, σ, u m , T only, sup t∈[0,t] E[(δX t ) 2 ] ≤ CE[X 2 0 ] δu L 1 (]0,T [) The continuity of u → T 0 E[(X t -X d ) 2
]dt follows and existence is then shown as for the deterministic case.

Existence of solution: the dynamic case

To establish existence by probabilistic methods proved to be rather hard. On the other hand, the problem can be converted into an equivalent distributed control problem by using the Kolmogorov equation of the PDF of X t . We will show in section 7 that the problem has a solution when the expected value in the criteria and κ are multiplied by a cutoff function for large X and

U = {u ∈ L 2 W 1,1 (R + ) : u m ≤ u ≤ u M , ∂ t u ∈ L 2 (Q + ), lim Y →∞ | log X|u(X, t) ≤ ∞.} 2.4. Discretization Let h = T
M and t m = mh. For any f , f m denotes an approximation of f (mh). Define (Euler scheme)

X m+1 = X m (1 + h(r -κX m -u m ) + σδW m ) ( 10
)
where δW m = W (m+1)h -W mh . Note that positivity may not be preserved by this scheme, but we may use [START_REF] Kamalapurkar | Reinforcement Learning for Optimal Feedback Control[END_REF]:

X m+1 = X m exp (r - σ 2 2 )h -(κX m + u m )h + σδW m .
Alternatively, to avoid computing exponentials, we may consider the Milstein scheme:

X m+1 = X m 1 + h(r -κX m -u m - σ 2 2 ) + σδW m + σ 2 2 (δW m ) 2 . ( 11
)
Proposition 3. Let X denote the solution of the continuous Milstein scheme given by

Xt = X m 1 + (t -mh) r -κX m -u m -σ 2 2 + σ(W t -W m ) + σ 2 2 (W t -W m ) 2 .
Assume X0 = X 0 > 0. Then on the event

A h = X m < 1 2κh + r-κX m -u M κ , m = 0, . . . , M -1 the ( Xt ) t∈[0,T ] almost surely has positive paths.
Proof. We proceed by induction. Starting from 11, we want the discriminant of the second degree polynomial equation in W t -W m to be negative, namely after obvious simplifications (since σX m = 0)

1 -2 1 + r -κX m -u m -σ 2 2 > 0 or, equivalently X m < 1 2κh + r-κX m -u M κ .
It iso clear that P(A h ) ↑ 1 as h → 0 (very fast indeed). For more details on the positivity preserving property of the Milstein scheme, see [START_REF] Pagès | Numerical Probability: An Introduction with Applications to Finance[END_REF]Chapter7.5].

Brute force solution of the non-dynamic control problem by Monte-Carlo

Recall that numerically δW m = √ hN m 0,1 where N m 0,1 denotes the numerical realization at time step m of a normal Gaussian random variable. For each of the 3 previous schemes there is an ζ h such that

X m+1 = ζ h (X m , u m , N 0,1 ), X 0 given, (12)
Let {X m k } K 1 be K realizations of [START_REF] Chollet | Deep learning with Python[END_REF]; by the law of large numbers :

E[(X t -X d ) 2 ] |mh ≈ E K [(X m -X d ) 2 ] := 1 K K 1 (X m k -X d ) 2 .
Problem (1) is approximated by min

um≤u m ≤u M J K := h M -2 0 E K [(X m+1 -X d ) 2 ] -αu m + β h 2 |u m+1 -u m | 2 ]
, with [START_REF] Chollet | Deep learning with Python[END_REF]. [START_REF] Bris | Existence and uniqueness of solutions to Fokker-Planck type equations with irregular coefficients[END_REF] Note that {u m } need not be Markovian; it may depend also on X p , p > m.

Solution with a gradient method

Gradient methods update {u m } M -1 1 by u m ← u m -µG m , where µ, the step size, is an appropriate scalar and G m = ∂ u m J K . Let X m k be computed by [START_REF] Brochier | Can overexploited fisheries recover by self-organization? Reallocation of the fishing effort as an emergent form of governance[END_REF]; let us differentiate this equation and multiply it by an adjoint variable P m k and sum:

M -2 m=0 -P m k δX m+1 k + P m k (δX m k (1 + h(r -κX m k -u m ) + σδW m k ) -X m k h(κδX m k + δu m )) = 0.
As δX 0 = 0, if we let P M -1 = P M -2 = 0, and define

v m k = 1 + h(r -2κX m k -u m ) + σdW m k , then - M -2 0 hP m k X m k δu m = M -2 1 (P m-1 k δX m k -P m k δX m v m k ) = M -2 1 δX m k (P m-1 k -P m k v m k ).( 14
)
Let us define P m k by

P m-1 k = P m k v m k -2(X m k -X d ), m = M -2, M -3, .., 1. Then M -1 m=1 δE K [(X m -X d ) 2 ] = 2 K K k=1 M -1 m=1 (X m k -X d )δX m k = - 1 K K k=1 M -1 m=1 (P m-1 k -P m k v m k )δX m k = 1 K K k=1 M -2 m=0 P m k X m k δu m k . ( 15
)
So the discrete gradient is

G m = 1 K K 1 P m k X m k -α -2 β h 2 (u m+1 -2u m + u m-1 ).

Implementation and examples

We use the L-BFGS-B function of the Python library minimize. It is an accelerated projected gradient method for problems with box constraints. We took 50 time steps and chose to solve a deterministic case with

σ = 0, K = 1, X d = 1, κ = 1.2, X 0 = 0.9, c = 1, a = 1, u m = 0.5, u M = 1. ( 16
)
The computer program is in file fishingQuotaBruteforce.py. The results with α = β = 0 are shown on the left in figure 1 with a saturation of the constraints u = u m when m < 7; but the solution is not bang-bang everywhere and u = 0.8 gives X m+1 = X m = X d , P m = 0. On the right in figure 1 the results are shown for α = 0.01, β = 0.001. On figure 2 a stochastic case of ( 13) is solved with the same parameters ( 16) except σ = 0.

Solution of the non-dynamic control problem by supervised learning

To reduce the CPU cost of solving

min u∈U J(u) := T 0 [ 1 K K 0 (X k (t) -X d ) 2 -αu(t) + β(∂ t u(t)) 2 ]dt (17) with dX k = X k (t)((r -κX k (t) -u(t))dt + σdW k (t), X(0) = X 0 , (18)
we shall train a neural network with input parameters [X 0 , r, κ] and output an optimal {u m } M -1

1
for these parameters, computed by the brute force method above. There is no need to vary the parameter X d because we can always return to the case X d = 1 by scaling. To build a training set we use a large number of synthetic solutions obtained with the brute force method of the previous section. However as noise is known not to affect neural network we use a low value K = 100 to generate 10000 solutions (it takes about 90 minutes on an intel core i7 machine). It works as follows:

• Choose [X 0 , r, κ] uniformly random in [X 0 m , X 0 M ] × [r m , r M ] × [κ m , κ M ].
• Generate K • M normal gaussian random values for dW t .

• Compute a discrete optimal solution {u m } M -1 1 of (17) with K=100 by the L-BFGS-B Python module.

• Then train a neural network [X 0 , r, κ] → {u m } M -1 1 with these 10000 synthetic solutions.

• Finally check the precision by using the network's answer {u m } M -1

1
to a new set of values for [X 0 , r, κ] and compare with the brute force solution corresponding to these parameters .

For an introduction to Statistical Learning with Neural Networks see [START_REF] Goodfellow | Deep Learning[END_REF]. We have used Keras (see [START_REF] Chollet | Deep learning with Python[END_REF]) to generate a 4-layer Neural Network (see figure 3). The input layer has 3 inputs; the two hidden layers have 50 neurons each and the output layer has M nodes. All layers use the ReLU activation. The two hidden layers are separated by a batch-normalisation. Notice that the error is about the same at all time steps. In figure 4 the performance of the trained Neural Network performance on 6 new cases are shown. Each case corresponds to a set of parameters [r, κ, X 0 ] and a new realization of the random noise. The "exact" control u is computed by the brute force method with K = 100 and compared with the solution u N N given by the neural network (values at the 50 time steps). In all cases except case 4 the control u N N given by the network is close to the "exact" one. But even in case 4, X N N computed by [START_REF] Brochier | Can overexploited fisheries recover by self-organization? Reallocation of the fishing effort as an emergent form of governance[END_REF] with u N N is close to X computed with u "exact".

Bellman's Stochastic Dynamic Programming for the dynamic problem

Consider the value function

V (t, X) = min u∈U T t E (X τ -X d ) 2 -αu + β|∂ t u| 2 dτ : dX τ = X τ dτ (r -κX τ -u τ ) + σX τ dW τ , X(t) = X} (19) Let [0, T ] = ∪ M 1 [(m-1)h, mh]. Let ζ h (X, u m , z
) denote one iterate of a numerical scheme [START_REF] Brochier | Can overexploited fisheries recover by self-organization? Reallocation of the fishing effort as an emergent form of governance[END_REF] [START_REF] Yong | Stochastic Controls Hamiltonian Systems and HJB Equations Application of Mathematics series[END_REF] for the SDE starting at X, for instance with [START_REF] Brochier | Can overexploited fisheries recover by self-organization? Reallocation of the fishing effort as an emergent form of governance[END_REF]

ζ h (X, u, z) = X + X(r -κX -u)h + Xσz √ h, z being the realization of a N 0,1 r.v. Let v M (X) = 0 and v m (X) = min u∈U (m+1)h mh E |ζ h (X, u, z) -X d | 2 -αu + β|∂ t u| 2 dτ + E[v m+1 (ζ h (X, u, z)] ≈ min u∈U hE |ζ h (X, u, z) -X d | 2 -αhu + β h |u m+1 (X) -u| 2 + E[v m+1 (ζ h (X, u, z)] . (20) Evidently E |ζ h (X, u, z) -X d | 2 = [X -X d + hX(r -κX -u)| 2 + hX 2 σ 2 .
To compute E[v m+1 (ζ h (X; u; z))], we use a quadrature formula with Q points {z q } Q 1 and weights {w q } Q 1 based on optimal quantization of the normal distribution N (0, 1) (see [START_REF] Pagès | Optimal quadratic quantization for numerics: the Gaussian case[END_REF] or [4, Chapter 5] and the website www.quantize.maths-fi.com for download of grids) so that

E[v m+1 (η h (X; u; z))] ≈ Q q=1 w q v m+1 (ζ h (X, u, z q )
Finally at every time step and every X j = j L/J, j = 1..J, with L >> 1, the result is minimized with respect to u ∈ U by a dichotomy. In this fashion {u m (X j )} J j=1 is obtained and a piecewise linear interpolation is constructed to prepare for the next time step u m-1 . With the same parameters as above except β = 0.1, L = 3, J = 50 and M = 50, the results of figures 5, 6 have been obtained. Once {u m (X)} M -1 0 is known, for any realization of the fishing model we use it as follows:

X m+1 = X m + hX m (r -κX m -u m (X m )) + X m σ √ hN m 0,1 , X 0 = X 0 .
Two such simulations are shown on figure 7,8. The results are not as good as before, meaning that the dynamic setting gives a poor solution to the control of the fishing site. These results should be compared with those of Appendix B which deals with a very similar problem for which there is a polynomial in log X solution when α = 0.

Solution with the Hamilton-Jacobi-Bellman partial differential equations

The Dynamic Programming equations can be written as a deterministic set of PDE: the Hamilton-Jacobi-Bellman (HJB) equations.

To this end we consider the same value function defined in (19)

Proposition 4. The value function verifies V (•, T ) = 0 and

0 = inf u∈U {(X -X d ) 2 -αu + β|∂ t u| 2 + ∂ t V + (r -κX -u)X∂ X V + σ 2 X 2 2 ∂ XX V }. ( 21
)

Proof

The proof given below is heuristic and non standard but it gives a link with the previous section. For the standard argument see for instance [START_REF] Yong | Stochastic Controls Hamiltonian Systems and HJB Equations Application of Mathematics series[END_REF]. We return to (20) above:

v m (X) = min u∈U (m+1)h mh E |ζ h (X, u, z) -X d | 2 -αu + β|∂ t u| 2 dτ + E[v m+1 (ζ h (X, u, z)] ≈ min u∈U E |X -X d | 2 -αu + β h 2 |u m+1 (X) -u| 2 h + E[v m+1 (ζ h (X, u, z)] (22)
As u M plays no role in the discretization we assume that u M = u M -1 . Hence when m = M -1 the minimization reduces to

v M -1 = h min u∈U {E |X -X d | 2 -αu} = h(|X -X d | 2 -αu m )
and the minimizer is

u M -1 = u m . When m < M -1 a similar argument gives, v m (X) ≈ min u∈U |X -X d | 2 -αu + β h 2 |u m+1 (X) -u| 2 h + E[v m+1 (ζ h (X, u, z)] . ( 23
)
With a Taylor expansion, let us approximate v m+1 (X + h) ≈ q M q=0 v (q) m+1 (X) h q q! . For clarity denote λ = r -κX -u, so

E[v m+1 (ζ h (X, u, z))] = v m+1 (X) + E[ q M 0 v (q) m+1 )(X)(X √ h) q (σz + λ √ h) q ] ≈ v m+1 (X) +Xv m+1 λh + X 2 2 v m+1 (λ 2 h + σ 2 )h + X 3 6 λσ 2 v (3) m+1 h 2 + X 4 24 σ 4 v (4) m+1 h 2 + o(h 2 ) E[v m+1 (ζ h (X, u, z))] ≈ v m+1 (X) + Xv m+1 (r -κX -u)h + X 2 2 v m+1 σ 2 h + o(h)
Let us use it in (23). Then u must solve

min -αu + β h 2 (u 2 -2uu m+1 ) -uXv m+1 , ⇒ u = u m+1 + h 2 2β (α + Xv m+1 ) (24)
but comply also with the box constraints, i.e.

u m = min u M , max u m , u m+1 + h 2 2β (α + Xv m+1 ) , and
v m ≈ v m+1 + h Xv m+1 (r -κX -u m ) + X 2 σ 2 2 v m+1 + |X -X d | 2 -αu m + β h 2 |u m+1 -u m | 2 . ( 25
)
As this corresponds to an explicit discretization of a parabolic PDE, h will have to be extremely small: h = O(δX 2 , where δX is the spatial discretization increment.

Numerical Results

We have implemented the method using FreeFEM++ [START_REF] Hecht | New development in FreeFem++[END_REF] with a slight change in (25): the term

X 2 σ 2 2 v m+1 is changed to X 2 σ 2 2 v
m so as to increase stability.

Numerical results shown on figure 9, 10 are obtained with the program in file pagesbellmanstoch.edp and the following parameters:

r = 2, κ = 1.2, σ = 0.3, X d = 1, α = 0.01, β = 0.1, u m = 0.5, u M = 1.
The computational domain is X ∈ (0, 3) divided in 400 intervals, t ∈ (0, 2), discretized with 40 time steps. Once {u m (X)} M -1 0 is known, for any realization of the fishing model we use it as follows:

X m+1 = X m + hX m (r -κX m -u m (X m )) + X m σ √ hN m 0,1 , X 0 = X 0 .
Two such simulations are shown on figure 11,20. Proof: With U = {u ∈ [u m , u M ]} and u m > 0, X being positive, the optimal control u * must be such that if there are 0 < t 1 , t 2 < T such that for all t ∈ (t 1 , t 2 ), u(X, t) is constant then 1. u * = u m when X∂ X V + α < 0, leading to ∂ t u * = 0 and

0 0.5 1 1.5 2 0.5 1 1.5 time u(X t , t), X t X t with quota quota u t X t without quota
∂ t V + (r -κX -u m )X∂ X V + σ 2 X 2 2 ∂ XX V = αu m -(X -X d ) 2 2. u * = u M when X∂ X V + α > 0 , leading to ∂ t u * = 0 and ∂ t V + (r -κX -u M )X∂ X V + σ 2 X 2 2 ∂ XX V = αu M -(X -X d ) 2 3. Otherwise when u m < u * (X, t) < u M , we must have V (T ) = 0 and 0 = (X -X d ) 2 -αu * + ∂ t V + (r -κX)X∂ X V + σ 2 X 2 2 ∂ XX V, X∂ X V + α = 0
The third option being impossible, the solution of the problem is bang-bang.

Remark 1. Note that this HJB dynamic solution is different from the non-dynamic solution of figure 1 6.3. Discretization One way to find u * when it is neither u m nor u M is to discretize in time and solve

inf u m g(u m ) := {-αu m + β h 2 |u m -u m+1 | 2 -u m X∂ X V }
The stationary point is at

-α + 2β h 2 (u m-1 -u m ) -X∂ X V m-1 = 0, i.e. ūm-1 = u m + h 2 2β (α + X∂ X V m-1
).

So we have to compare g(u m ), g(ū m ) and g(u M ) and choose the smallest. The end result is very similar to the previous section. Two simulations are done with β = 0 and shown on figure 13,14. The control is indeed bang-bang. 0 0.5 1 1.5 2 0.5

X t with quota quota u t X t without quota

X t with quota quota u t X t without quota

Solution with the Kolmogorov equation

The Kolmogorov equation for the probability density ρ(X, t) of {X t } T 0 is:

∂ t ρ + ∂ X [(r -κX -u)Xρ] -∂ XX [ X 2 σ 2 2 ρ] = 0, ρ(X, 0) = ρ 0 (X), ∀X ∈ R + . (26)
The optimal control problem (1) is:

min u∈U J(u) := Q + (X -X d ) 2 -αu(X, t) + β|∂ t u(X, t)| 2 ρ(X, t)dXdt : with (26). ( 27
)
Integration is done on Q + := R + ×]0, T [. The conditions for having equivalence between the stochastic control written as ( 1) and ( 27) are detailed in [START_REF] Bris | Existence and uniqueness of solutions to Fokker-Planck type equations with irregular coefficients[END_REF]. To avoid using weighted Sobolev spaces we make a change variable y = log X, set v = r -σ 2 2 -κe y -u(e y , t) so that dY t = vdt + σdW t and then study the existence of a PDF for Y t given by

∂ t ρ + ∂ Y (vρ) -∂ Y Y [ σ 2 2 ρ] = 0, ρ(Y, 0) = ρ 0 (Y ), ∀Y ∈ R, ∀t ∈]0, T [. ( 28
)
To analyze the last term in the criteria we write it as

β Q R |∂ t v n | 2 ρ n = β Q R |∂ t v n | 2 ρ + β Q R |∂ t v n | 2 (ρ n -ρ).
By the lower semi-continuity of semi-norms,

lim n→∞ Q R |∂ t v n | 2 ρ ≥ Q R |∂ t v| 2 ρ. Finally the last term above is bounded by β (∂ t v n ) 2 L 2 (Q R ) ρ n -ρ L 2 (Q R ) ≤ βK ρ n -ρ L 2 (Q R ) → 0. Consequently, the criteria in (30) is l.s.c. Remark 2. If + R X 2 ρ 0 (X)dX < ∞, it
may be possible to let R → ∞ in the proof above.

Computation of gradients

Consider the variational form of the Kolmogorov equation: find ρ ∈ L 2 (0, T, V ) such that, for all ρ ∈ V ,

R + ρ∂ t ρ -(r -κX -u)Xρ∂ X ρ + σ 2 2 ∂ X (X 2 ρ)∂ X ρ = 0, ρ(0) given. ( 32
)
Calculus of variations leads to δρ(0) = 0 and

R + ρ∂ t δρ -(r -κX -u)Xδρ∂ X ρ + σ 2 2 ∂ X (X 2 δρ)∂ X ρ = - R + ρX∂ X ρδu. (33)
Define the adjoint p by p (T ) = 0 and, for all p ∈ V :

R + p∂ t p + p(r -κX -u)X∂ X p - σ 2 2 ∂ X (X 2 p)∂ X p +p1 Ω [(X -X d ) 2 -αu + β|∂ t u| 2 ] = 0. (34)
Add ( 33) and (34) with ρ = p and p = δρ. It gives

R + ∂ t (p δρ) + Ω [(X -X d ) 2 -αu + β|∂ t u| 2 ]δρ = R + ρX∂ X p δu.
As p (T ) = 0 and δρ(0) = 0, an integration in time gives

Ω×]0,T [ [(X -X d ) 2 -αu + β|∂ t u| 2 ]δρ = - Q + ρX∂ X p δu.
Finally, from (27),

δJ = < grad u J, δu >= Ω×[0,T ] (X -X d ) 2 -αu + β|∂ t u| 2 δρ -αρδu + 2ρ∂ t u∂ t δu . = - Q + ρX∂ X p δu + Ω×]0,T [ ρ[-αδu + 2∂ t u∂ t δu] (35)

Results

A steepest descent algorithm with fixed step size and projected gradient on the box constraints is

v = u k -µgrad u J k , u k+1 = max{u m , min{u M , v}}. (36)
Note that we have assumed that ∂ t u(0) = ∂ t u(T ) = 0. The same parameters as above have been used except that the computational domain is (0, 5) instead of (0, 3). The PDF of X 0 is a Gaussian curve with σ 0 = 0.2. The computer program is in the file GPOP2dxt.edp. After 100 iterations J is a fifth of its initial value. Figure 15 shows the surfaces u,P and ρ versus X, t. Note that (X, t) → u(X, t) is either 0.5 or 1, implying that all controls u(X t , t) will be bang-bang. Figures 16 and17 show two applications, one with X 0 = 0.8 and the other with X 0 = 1.2. The method achieves a reasonable stabilization of the fish biomass, but the quotas vary rapidly, as if the β -term wasn't playing its regularizing role. Note that (X, t) → u(X, t) is either 0.5 or 1, implying that all controls u(X t , t) will be bang-bang. There is a vertical scale factor of 0.5 applied.

Solution by Itô calculus

In computational finance Itô calculus is preferred over the Kolmogorov setting, primarily because the initial condition for the PDF could be singular. Consider again Problem [START_REF] Bellman | Dynamic Programming[END_REF]. Itô calculus applied to

j t = T t E (X τ -X d ) 2 + β(∂ t u) 2 -αu dτ
as a function of X t gives:

d t j = ∂ t jdt + E[∂ X jdX] + 1 2 E[∂ XX jdXdX] = dt ∂ t j + ∂ X j(r -κX -u)X + X 2 σ 2 2 ∂ XX j = -dt[(X -X d ) 2 + β(∂ t u) 2 -αu] (37)
This formulation is dual to the Kolmogorov calculus presented above. Let us compute j as the solution at τ = 0 and X = X 0 of the PDE, backward in time and initiated by j(X, T ) = 0, ∀X ∈ R + , τ ∈]0, T [:

∂ t j + (r -κX -u)X∂ X j + X 2 σ 2 2 ∂ XX j = -[(X -X d ) 2 + β(∂ t u) 2 -αu] (38)
The optimal control problem ( 1) is now min u∈U R j(X, 0)ρ 0 (X)dX : subject to (38)

where ρ 0 (X) is the PDF of X 0 . We have to restrict X < X M for the right hand side of the PDE to be in L 2 (R).

On Ω := (0, X M ) the variational formulation of the PDE consists in finding j in a weighted Sobolev space:

j ∈ V := {g ∈ L 2 (Ω) : X∂ X g ∈ L 2 (Ω)}, j(•, T ) = 0, Ω (∂ t j + (r -κX -u)X∂ X j) ĵ - Ω σ 2 2 ∂ X (X 2 ĵ)∂ X j + σ 2 2 X M j(X M ) ĵ(X M ) = - Ω [(X -X d ) 2 + β(∂ t u) 2 -αu] ĵ, ∀ ĵ ∈ V. (40)

Gradient computation

The gradient of j with respect to u results from the following calculus:

Ω×]0,T [ grad u jδu := Ω ρ 0 δj(X 0 , 0) subject to Ω (∂ t δj + (r -κX -u)X∂ X δj) ĵ - Ω σ 2 2 ∂ X (X 2 ĵ)∂ X δj + σ 2 2 X M δj(X M ) ĵ(X M ) = Ω ĵ[(X∂ X j + α)δu -2β∂ t u∂ t δu], ∀ ĵ ∈ V. (41)
Define p by

Ω (p∂ t p -(p(r -κX -u)X∂ X p) + Ω σ 2 2 ∂ X (X 2 p)∂ X p - σ 2 2 X M p(X M )p(X M ) = 0, p(0) = ρ 0 (42)
Then adding (42) to (41) with p = δj and ĵ = p and remembering that δj(X, T ) = 0, leads to

Ω (p∂ t δj + δj∂ t p) = Ω×]0,T [ [p(X∂ X j + α)δu -2β∂ t u∂ t δu] (43)
On the other hand

- Ω ρ 0 δj(X, 0) = Ω δj(X, T )p(T ) - Ω δj(X, 0)p(0) = Ω×]0,T [ ∂ t (pδj) (44) Consequently Ω×]0,T [ grad u jδu = Ω ρ 0 δj(X 0 , 0) = - Ω×]0,T [ [p(X∂ X j + α)δu -2β∂ t u∂ t δu]
Hence this method is a minor variation of the method with the Kolmogorov equation: ρ and p have similar equations and p and j too.

Limit with vanishing volatility

We have seen earlier that the optimal solution of the non-dynamic control is very close to the solution of the deterministic control problem.

We explore here the convergence of the PDE method to the deterministic case by letting σ tend to zero. For clarity we take α = β = 0. Consider again (39). Let Y (X, τ ) be the solution of

Ẏ (τ ) = (r -κY (τ ) -u(Y (τ ), τ )) Y (τ ), Y (t) = X; (45)
then

∂ t j + (r -κX -u)X∂ X j = d dτ ρ(Y (τ, X), τ )| τ =t . (46)
Hence ( 38) is rewritten as

d dτ j(Y (τ, X), τ )| τ =t = - σ 2 2 ∂ XX (X 2 ρ) -(X -X d ) 2 .
Consequently, when σ → 0, j becomes equal to j, the integral of -(X -X d ) 2 on the streamline (45) which passes through X 0 at t = 0:

j(Y (T ), T ) -j(Y (0), 0) = -j(X 0 , 0) = - Y (T ) Y (0) (X(s) -X d ) 2 i.e. j(0) = T 0 (X(t) -X d ) 2 dt : Ẋ = (r -κX -u)X, X(0) = X 0 . (47)
This analysis shows that HJB becomes numerically singular when σ is too small.

Concluding discussion

Let us compare the 3 classes of methods:

1. Non-dynamic control (NDC) solved by Monte-Carlo. 2. Dynamic Programming in a stochasting setting (DP).

Dynamic Programming solved with Partial Differential equations (HJB).

NDC is a mathematical oddity because it does not take into account the necessity of a Markovian solution forbidding an attempt to use the current noise level to predict the next step (feedback); but from the practical point of view it is certainly the best solution for the fishing site regulator, especially when it is coupled with Supervised Learning. After all, provided that the model fits the reality, a Markovian solution is not needed. By reinitializing the data every so often at {T k } J 0 and solving J control problem on smaller time intervals [T j-1 , T j ], on could borrow from Model Predictive Control (see [START_REF] Boyd | Linear Controller Design -Limits of Performance[END_REF]). Yet, it can be seen on figure 4 that the precision is hardly better at early times than at later times. On this fishing site control problem DP and HJB are different implementations of the same approach. All four implementations gave disappointing results, probably because the parameter chosen gave a solution close to a bang-bang saturation of constraints and not far from the elementary but intuitive strategy which consists in taking u t = u m when X t > X d and u M otherwise. Another possible reason is that penalization of |∂ t u| works for NDC but doesn't for DP and HJB. A penalization of |∂ X u| is probably needed as well. HJB with PDEs is the most taxing method numerically and DP is much faster and also capable of generating approximate solution which do not saturate the constraints as much. NDC requires brute force Monte-Carlo, which is also very expensive, yet with neural network training for NDC outperform traditional HJB and DP methods. This problem belongs to a class for which the user cannot measure precisely the random errors in the data and the model; the number of boats at sea and the quantity of fish caught everyday may just reflect the roughness of the model rather than the stochasticity of the measurements; nevertheless this point needs to be investigated further and in the mean time the non-dynamic control problem makes sense and gives better results. Furthermore, when coupled with a neural network the optimization process is very fast, once the learning phase is done. Finally, porting the trained network to a hand held device is also a major advantage.

Reformulation

To prevent fish extinction a constraint is set on the total catch min(qX(t), Q(t))E(t)). Optimization will then define Q(t). It is expected that Q(t) < qX(t) otherwise the policy of quota is irrelevant in the sense that the fisherman is given a maximum allowed catch which is greater than what he could possibly catch. So let us add a constraint Q(t) ≤ Q M so that min(qX(t), Q(t)) = Q(t). Denote u(t) = E(t)Q(t)/X(t). As said before it is not feasible too impose to small a quota, so finally U = {u m ≤ u(t) ≤ u M }. where u M is less than qE(t) at all time.

12.1. More constraints on quotas Small quotas are unpopular. To avoid them we may use penalty and add to the criteria -T 0 αu(t). More over, it will be unpopular if u(t) = u m one day and u(t) = u M the next day, so to avoid too many changes we may add to the criteria a term like β| du dt | 2 . An optimal policy is a solution of min Q∈Q J(u) :=

Figure 1 shows the solution with the same parameters as for figure 18. There are differences because it is a different set of constraints. The policy given to each fisherman is Q(t) for the full system and u(t)X(t)/E(t) for the reduced system.

Proposition 1 .

 1 Problem (5) has a solution.

Figure 1 :

 1 Figure 1: Solution of (1) discretized by[START_REF] Bris | Existence and uniqueness of solutions to Fokker-Planck type equations with irregular coefficients[END_REF], versus m, in the deterministic case (σ = 0, K = 1) with β = 0 (left), and β = 0.001 (right).The blue curve is the control u; the biomass X t is in ochre and the fishing effort in green.

  3. The results are shown for K = 1, 10, 100, 10000 Monte-Carlo samples. Note that the brute force solution for K >> 1 is quite close to the noiseless solution of figure 1-b.

Figure 2 :

 2 Figure 2: Solution of (13) in the stochastic case with σ = 0.3: convergence of the Monte-Carlo approximation . Left to right, up first: K=1, K=10, K=100,K=1000, K=10000. Already the results are fairly converged at K=100. In the last figure the optimal control is computed with K=100 and used to solve[START_REF] Balakrishnan | Adaptive-critic-based neural networks for aircraft optimal control[END_REF] with a new random noise, thereby showing that the solution is still acceptable.

Figure 3 :

 3 Figure 3: The Neural Network has 3 inputs, two hidden layers each with 50 neurons (10 shown here) and separated by batch-normalisation and an output layer with M output (M=15 shown here) . Right: typical convergence curve of the loss function during the training phase.

4. 1 . 3 )

 13 A numerical testWe have taken 50 time steps and K=100 for the brute force part. The following parameters are used:T = 2, α = 0.01, β = 0.001, σ = 0.3, X d = 1, u m = 0.5, u M = 1For the Neural Network training we have used 10000 samples in batches of 32 and a maximum of 2000 epochs. The input parameters are r = np.random.uniform(1.9,2.2) kappa = np.random.uniform(1,1.3) X0=np.random.uniform(0.7,1.Finally for the fishing effort equation, a = c = 1. The computer program is written in Python and is in the file fishingQuotaIA4.py. The learning phase produced the following precision on each of the time steps :Test loss: 0.03751661651 Average relative error |1-u_NN/u| for the 49 time steps = [0.03185646 0.03954296 0.04612997 0.03440035 0.03665108 0.04323719 0.03549499 0.04364510 0.04247633 0.03737225 0.03448902 0.03950258 0.04441193 0.04149515 0.03366303 0.04054435 0.03978632 0.04160346 0.0404573 0.03319214 0.03800414 0.03229879 0.04076255 0.03475258 0.03333832 0.03761944 0.04427279 0.03837532 0.03909892 0.03384874 0.03707403 0.03838878

Figure 4 :

 4 Figure 4: On each of the 6 plots X and u are shown: True discrete solution in green for X and red for u. NN solution in blue for X N N and ochre for u N N . The plots correspond to the following values of the parameters: r = 2.0228569911704177, κ = 1.0203759108766386, X 0 = 1.1122268356898797, r = 2.0688190264082973, κ = 1.2641632371604865, X 0 = 0.7852264857858362, r = 1.9151933649036637, κ = 1.0141074063562392, X 0 = 0.8299957532445648, r = 2.0738603520220655, κ = 1.187620910330912, X 0 = 1.1344043926097065, r = 2.0463252714459936, κ = 1.122409869546898, X 0 = 0.7116598707629698, r = 1.9002636936181148, κ = 1.189933237216782, X 0 = 0.7550552945829774.

Figure 5 :

 5 Figure 5: Solution of the Dynamic Programming equation: u(x, t).

Figure 6 :

 6 Figure 6: Solution of the Dynamic Programming equation: V (x, t).

Figure 7 :

 7 Figure 7: Simulation of the fishing model with a quota function computed by Dynamic Programming and X 0 = 0.8.

Figure 8 :

 8 Figure 8: Simulation of the fishing model with a quota function computed by Dynamic Programming and X 0 = 1.2.

Figure 9 :Figure 10 :

 910 Figure 9: Solution of the Dynamic Programming equation: side view.

Figure 11 :

 11 Figure 11: Simulation of the fishing model with a quota function computed by Dynamic Programming and X 0 = 0.8.

Figure 12 :

 12 Figure 12: Simulation of the fishing model with a quota function computed by Dynamic Programming and X 0 = 1.2.

6. 2 .

 2 Analysis when β = 0 Proposition 5. When β = 0, the solution is always bang-bang.

Figure 13 :

 13 Figure 13: Simulation of the fishing model with a quota function computed by HJB when β = 0 and X 0 = 0.8.

Figure 14 :

 14 Figure 14: Simulation of the fishing model with a quota function computed by HJB when β = 0 and X 0 = 1.2.

Figure 15 :

 15 Figure 15: From left to right: optimal ρ,p ,u. Horizontal axis is X∈]0, X M [, vertical is t ∈ [0, T ] with X M = 5, T = 2, X d = 1, u m = 0.5, u M = 1.Note that (X, t) → u(X, t) is either 0.5 or 1, implying that all controls u(X t , t) will be bang-bang. There is a vertical scale factor of 0.5 applied.

Figure 16 :

 16 Figure 16: Optimal biomass and quota function computed by Optimal Control and compared with the solution with u = u M when X 0 = 0.8.

Figure 17 :

 17 Figure 17: Optimal biomass and quota function computed by Optimal Control and compared with the solution with u = u M when X 0 = 1.2.

Figure 18 :

 18 Figure 18: Left: Quota (blue) and fishing effort (ochre). No quota means that Q(t) = 1 for all t. Right: biomass with (blue) and without quota (ochre).

  Then we need u(t) ≤ qE(t) and the problem becomesmin t) -X d (t)) 2 dt : dX dt = X (r -κ X -u(t)) , X(0) = X 0 } (53)In the paper, the optimization is performed without computing E; then E is computed with the optimal u by solving dE dt = a -cE -uX, E(0) = E 0 .

T 0 |X -X d | 2 -

 02 αu + β| du dt | 2 : subject to dX dt = X (r -κ X -u(t)) , X(0) = X 0 .

Acknowledgements

Our sincere gratitude goes to François Murat for his help to minimize the hypotheses necessary for the existence of a solution to the Kolmogorov distributed control problem. The computer programs in Python and FreeFEM++ can be obtained by email from olivier.pironneau@sorbonne-university.fr.

The PDE is understood in the sense that for all ρ ∈ H 1 (R),

Note that v can be changed to v -L, L ∈ R by changing ρ to ρe -Lt . Hence when lim Y →∞ u = u ∞ = 0 we can take L = r -σ 2 2 -u ∞ and assume that lim Y →∞ v = 0.

, given by (28).

Consequently, for the problem to make sense, (r -κX -u) needs to be multiplied by a regular cutoff function, tending to zero fast enough when X → +∞. We shall henceforth assume that κ is a function of X which is zero when X > X M , for some X M .

and consider

Problem (30) has a solution.

Proof: Note first that J ≥ αv m . Note also that when u M is large enough so that

is finite and V is non empty the infimum is finite. Multiplying the PDE by ρ and integrating over R leads to

From the maximum principle we also know that 0 ≤ ρ n ≤ 1 in Q and that R ρ(Y, t)dY = 1 a.e. Now let us take a minimizing sequence {v n } n∈N for (30). Then

It follows from above that there is a v such that, for a subsequence, v n → v inL ∞ -weak-star and weakly in L 2 (Q) and a ρ such that

By the weak convergence of v n and the strong convergence of

APPENDIX A 11. A model with fishing quota

Let X(t) be the fish biomass at time t, E(t) the fishing effort -interpreted as the number of boats at sea -and Q(t) is the fishing capacity per boat. In [START_REF] Brochier | Can overexploited fisheries recover by self-organization? Reallocation of the fishing effort as an emergent form of governance[END_REF] Q(t) = qX(t),with the catchability q constant, meaning that the more fish there is the more fishermen will catch them. In this study a quota Q(t), given to each fisherman, is imposed on the maximum weight of fish caught on a day t ; hence the total amount of fish caught on a day is min(qX(t), Q(t))E(t). The modified model is an equation for X(t) which says that the biomass changes due to the natural growth/decay rate r, the long time limit κ of X and the depletion due to fishing:

Biomass dynamics. (48)

Let F (t) = min(qX(t), Q(t)). The fishing effort is driven by the profit pF minus the operating cost of a boat c, where p being the price of fish.

Fishing effort driven by profit,

The price is driven by the difference between the demand D(p) and the resource F E:

where Φ is the inverse time scale at which the fish market price adjusts. When Φ << 1, (50) may be approximated by

With a = a /γ and κ = κ γ, the whole system is :

Now notice that by changing the initial conditions X 0 → X 0 /γ we may take γ = 1.

Remark 3. Finally with t = tq, Q = Q/q and (r, κ", ã, c) = (r, κ, a, c)/q the above system is identical but now q = 1.

Remark 4. The ODE system is of the form

Global existence of solution on [0, T ] is not known because ψ is only locally Lipschitz and not sublinear for large E, X.

A numerical simulation with and without optimal quota is shown on figure 18 with the following parameters (see file fishingQuota0.py and for the definition of α and β, see below):

Notice that q is big enough so that min(qX(t), Q(t)) = Q(t), ∀t. Notice also that without quota the biomass tends to zero, in this case.

APPENDIX B

13. An analytical solution for a similar problem

It was said earlier that Y t = log X t could be used to simplify the problem. As the goal is to reach

If the constraints are not active, the equation for the optimality of u is

When α = 0 there is a polynomial solution to this system: V = aY 2 +bY +c, v = KY +k where a, b, c, K, k are functions of time with a(T

One must also check that u