Fany Sardenne

Nathalie Bodin

Leo Barret

Laura Blamey

Rodney Govinden

Kettyna Gabriel

Rosabella Mangroo

Jean-Marie Munaron

François Le Loc'h

Antoine Bideau

Fabienne Le Grand

Magali Sabino

Paco Bustamante

David Rowat

Diets of spiny lobsters from Mahé Island reef, Seychelles, inferred by trophic tracers

Keywords: Fatty acids, Stable isotopes, Bayesian mixing model, Panulirus sp., tropical reef habitats, Western Indian Ocean

Spiny lobsters (Panulirus longipes, P. penicillatus and P. versicolor) are an important resource in Seychelles, where they inhabit coastal carbonate and granite reefs that have been impacted by multiple coral bleaching events over the past two decades. Little is known about their biology and ecology in this region. Interspecific competition for food resources was previously suggested, but no quantitative data on the diet of spiny lobsters were available. Using carbon and nitrogen stable isotope compositions and fatty acid profiles of three spiny lobster species and their potential prey, a Bayesian mixing model for diet estimation was applied to compare the diet proportions of spiny lobsters among species and between reef types (carbonate and granite reefs). Model outputs suggested the three lobster species consume mainly crustaceans (Anomoura hermit crabs; half of the diet), then Echinoidea (sea urchins), algae and molluscs. P. versicolor was found to consume slightly more molluscs and algae than the two other studied species, which was consistent with its lower trophic level (2.4 vs 2.8 for the two other species). Trophic level did not increase with carapace length of spiny lobsters, but large individuals had higher carbon isotopic values suggesting that they might feed closer to the coast or more on detritus feeders than their smaller congeners. Diets of spiny lobsters were fairly similar between carbonate and granite reefs, except that lobster inhabiting granite reefs consumed more sea urchins. While our overall findings were consistent with gut contents of Panulirus spp. from other world regions, they should be confirmed, as the discrimination of several prey based on trophic tracers was low, which increased mixing model uncertainty.

Introduction

To maintain ecosystems in a healthy, productive, and resilient condition, several elements affecting species populations should be considered, such as interactions with other species or effects of environmental changes [START_REF] Pikitch | Ecosystem-Based Fishery Management[END_REF]. To reach this goal, ecosystem-based fisheries management requires fisheries research to account for relationships between target species and their habitats, including trophic interactions [START_REF] Bellchambers | Abundance and size of western rock lobster (Panulirus cygnus) as a function of benthic habitat: implications for ecosystem-based fisheries management[END_REF][START_REF] Tam | Towards ecosystem-based management: identifying operational food-web indicators for marine ecosystems[END_REF]. Spiny lobsters Panulirus longipes, P. penicillatus and P. versicolor are important ecological, social and economic resources in the Seychelles Archipelago, a small island state in the Western Indian Ocean (Fig. 1). They are found on shallow coral and granite reefs and in recent years, recurrent coral bleaching events have resulted in the loss of coral cover and structural complexity on these reefs and an increase in macroalgae cover (Graham et al., 2006;[START_REF] Harris | Scleractinian coral communities of the inner Seychelles 10 years after the 1998 mortality event[END_REF]. Although variable depending on the year, approximately three metric tons of spiny lobster (P. longipes, P. penicillatus and P. versicolor) are landed each year by the artisanal fishery (Barret, 2019), most of it being sold to local hotels and restaurants. Declines in catch per unit effort (CPUE) (24%) in the 2016-2017 season, and continued suppressed CPUE from independent fishery surveys have resulted in a temporary closure of the lobster fishery (SFA Fisheries Research Section, 2019), and highlighted the need for an ecosystem-based approach to fishery management, particularly given the new ecological conditions, i.e. recurrent coral bleaching events and likely habitat modification.

Despite the importance of spiny lobsters for the Seychelles, little is known about their biology and ecology in the region (Barret, 2019). As part of a precautionary approach to fisheries management, a number of management measures are in place, including the prohibition of traps to catch lobsters (they may only be caught using snorkelling) (Payet and Isidore, 2002). This avoids a possible modification of the diet toward bait, as has been observed in Western Australia [START_REF] Macarthur | Habitat surrounding patch reefs influences the diet and nutrition of the western rock lobster[END_REF]. However, the diet of spiny lobsters in the Western Indian Ocean remains unknown. A recent study based on trophic tracers, i.e. carbon and nitrogen stable isotope compositions (SI; δ 13 C and δ 15 N) and fatty acids (FAs) profiles, compared trophic niches of three spiny lobster species and found a high probability of dietary competition among these species (Sabino et al., 2020). While this study did not analyse the potential prey of spiny lobsters, results were consistent with the gut content analysis of Panulirus spp. in other world regions: Panulirus spp. are opportunists, mainly feeding on crustaceans and molluscs, and to a lesser extent on red algae [START_REF] Briones-Fourzán | Coral reef degradation differentially alters feeding ecology of co-occurring congeneric spiny lobsters[END_REF][START_REF] Castañeda-Fernández-De-Lara | Feeding ecology of juvenile spiny lobster, Panulirus interruptus, on the Pacific coast of Baja California Sur, Mexico. New Zeal[END_REF][START_REF] Joll | Natural diet and growth of juvenile western rock lobsters Panulirus cygnus George[END_REF][START_REF] Mashaii | Feeding habits of the scalloped spiny lobster, Panulirus homarus (Linnaeus, 1758) (Decapoda: Palinuridae) from the South East Coast of Iran[END_REF]. However, the proportion of algae might be underestimated due to a rapid gastric evacuation (<12 h; [START_REF] Waddington | Variation in evacuation rates of different foods skew estimates of diet in the western rock lobster Panulirus cygnus[END_REF]. Additionally, the predominant benthic habitat surrounding a reef can be an important source of variation for the diet of spiny lobsters [START_REF] Blamey | Diet of the spiny lobster Jasus paulensis from the Tristan da Cunha archipelago: Comparisons between islands, depths and lobster sizes[END_REF]Goñi et al., 2001;[START_REF] Macarthur | Habitat surrounding patch reefs influences the diet and nutrition of the western rock lobster[END_REF]. Similar variability in diet probably occurs for spiny lobsters from the Seychelles, as their potential for dietary competition was higher in granite reefs than in carbonate reefs (Sabino et al., 2020), probably in relation to the higher prey diversity in carbonate reefs (Robinson et al., 2019). However, quantitative diet studies are needed to better describe their trophic ecology while considering their different types of reef habitat in the Seychelles.

SI and FAs are intrinsic ecological tracers commonly used to infer assimilated diet in consumers [START_REF] Ramos | Trace me if you can: the use of intrinsic biogeochemical markers in marine top predators[END_REF], including in spiny lobsters (e.g. [START_REF] Blamey | Diet of the spiny lobster Jasus paulensis from the Tristan da Cunha archipelago: Comparisons between islands, depths and lobster sizes[END_REF][START_REF] Guest | Trophic effects of fishing southern rock lobster Jasus edwardsii shown by combined fatty acid and stable isotope analyses[END_REF]Sabino et al., 2020;Waddington et al., 2008). δ 13 C values allow for the discrimination of habitat types (e.g. coastal vs offshore, the latter being characterised by depleted δ 13 C values; [START_REF] France | Carbon-13 enrichment in benthic compared to planktonic algae: foodweb implications[END_REF] and δ 15 N values inform on the species' trophic level due to a 15 N enrichment with increasing trophic levels [START_REF] Vander Zanden | Comparing trophic position of freshwater fish calculated using stable nitrogen isotope ratios (δ 15 N) and literature dietary data[END_REF]. SI in muscle reflect the diet of spiny lobsters in the last three months [START_REF] Waddington | Diet quality and muscle tissue location influence consumerdiet discrimination in captive-reared rock lobsters (Panulirus cygnus)[END_REF]. FAs makeup lipids, and FAs of fat depots such as the hepatopancreas in crustaceans, reflect the FA profiles of prey assimilated by consumers in the last few days to weeks, depending on species [START_REF] Antonio | Influence of diet on the metabolic turnover of stable isotope ratios and fatty acids in the omnivorous shrimp Palaemon peringueyi[END_REF][START_REF] Shu-Chien | Effect of dietary lipid source on expression of lipid metabolism genes and tissue lipid profile in juvenile spiny lobster Sagmariasus verreauxi[END_REF]. Both SI and FA data can be used for the quantitative diet reconstruction through mixing models (e.g. [START_REF] Jankowska | Modification of benthic food web structure by recovering seagrass meadows, as revealed by trophic markers and mixing models[END_REF][START_REF] Neubauer | Bayesian estimation of predator diet composition from fatty acids and stable isotopes[END_REF][START_REF] Young | Investigating diet patterns of highly mobile marine predators using stomach contents, stable isotope, and fatty acid analyses[END_REF].

The objective of this study was to compare the trophic ecology of the three tropical spiny lobsters commonly caught in the Seychelles (P. longipes, P. penicillatus and P. versicolor) using SI and FA data from spiny lobsters and their potential prey. Specifically, we were interested in (1) what are the main prey of spiny lobster in Seychelles, (2) do diet and trophic level differ between spiny lobster species, and (3) does spiny lobster diet differ between reef types (carbonate vs granite reefs)?

Material and methods

Organism collection

Spiny lobsters and their potential prey were collected from the shallow coastal waters from reefs around the coasts of Mahé, the main island of the Seychelles Archipelago, home to around 90% of the country's population, and supporting 60% of the licensed lobster fishermen. Based on literature (Table S1), potential lobster prey (n=74 individuals from 13 species and 7 phylogenetic classes) were collected by scuba-divers at four sites (maximum 10 meters depth) (Fig. 1) in April 2018 and October 2019. Spiny lobsters (26 P. longipes, 31 P. penicillatus and 7 P. versicolor) were collected by divers on snorkel from 15 sites (12 granite and 3 carbonate reefs) in October 2018 and May 2019 (Table 1). Granite reefs are granitic rocky reefs with coral growth on the granite substrate, and carbonate reefs are continuous carbonate fringing reefs [START_REF] Harris | Scleractinian coral communities of the inner Seychelles 10 years after the 1998 mortality event[END_REF]. Note that many carbonate reefs have deteriorated or collapsed following multiple coral bleaching events (Obura et al., 2017) and therefore most lobsters in our study were found on granitic reefs. All organisms were frozen at -80 •C immediately after collection until further processing (less than seven days after collection) at the Seychelles Fishing Authority Research Laboratory. For spiny lobsters, a sample of hepatopancreas (for FA analysis) and a sample of tail muscle (for SI analysis) were carefully dissected from the organism. Potential prey were identified to species level whenever possible and their tissues were sampled for both FA and SI analysis. For echinoderms, a sample of muscle was carefully dissected from the Aristotle's lantern; for bivalves, gastropods and anomurans (decapod malacostracanshermit crabs), the whole-body tissue was collected after removing the shells; finally, polyps of Palythoa (anthozoan) were kept whole and algae (green algae of the Chlorophyta phylum; brown algae of the Ochrophyta phylum; and red algae or Rhodophyta) were kept with their small epiphytes (<1 mm).

All prey and spiny lobster samples were then kept at -80 °C for a maximum of five months, freeze-dried for 72 h, then ground to a homogeneous powder with a ball mill, and stored again at -80 °C prior to subsequent analysis.

Stable isotope analysis

For potential prey items, carbon isotope composition was analysed on a subsample acidified with HCl at 0.5 to 2 N (according to the effervescence) to remove carbonates, then rinsed three times with distilled water, and freeze-dried for 48 hours. Nitrogen isotope composition was analysed on a subsample of the bulk powder. For lobster muscle, the bulk powder was analysed for both carbon and nitrogen isotope composition, as the muscle did not contain carbonate neither high lipid content (fat content is 0.7% wet weight and C:N ratio remains < 3.5). Acidified and bulk powders were weighted into tin capsules and analysed by continuous flow on a Flash EA2000 elemental analyser coupled to a Delta V Plus isotope ratio mass spectrometer (Thermo Fisher scientific) at the Pôle Spectrométrie Océan, University of Brest, France, or on a Flash EA1112 elemental analyser coupled to a Delta V Advantage isotope ratio mass spectrometer (Thermo Fisher scientific) at the LIENSs stable isotope facility, University of La Rochelle, France. Calibrations were based on reference materials (IAEA-600, IAEA-CH-6 and IAEA-N-2). Results were reported in the δ unit notation and expressed as parts per thousand (‰) relative to the international standards (atmospheric N2 for nitrogen and Vienna-Pee Dee Belemnite for carbon). Analytical precision based on replicate measurements of acetanilide was < 0.15 ‰ for both δ 13 C and δ 15 N values at both laboratories. Inter-laboratory replicate measurements on reference materials (n=9), acetanilide (n=3), urea (n=3), casein (n=3), and lobster samples (n=3) found no significant differences in results between laboratories, with all differences below the analytical precision (< 0.15 ‰). C:N ratios were determined from % element weights.

Stable isotope analysis

For potential prey items, carbon isotope composition was analysed on a subsample acidified with HCl at 0.5 to 2 N (according to the effervescence) to remove carbonate, then rinsed three times with distilled water, and freeze-dried for 48 h. Nitrogen isotope composition was analysed on a subsample of the bulk powder. For lobster muscle, the bulk powder was analysed for both carbon and nitrogen isotope composition, as the muscle did not contain carbonate neither high lipid content (fat content is 0.7% wet weight and C:N ratio remains < 3.5). Acidified and bulk powders were weighted into tin capsules and analysed by continuous flow on a Flash EA2000 elemental analyser coupled to a Delta V Plus isotope ratio mass spectrometer (Thermo Fisher scientific) at the Pôle Spectrométrie Océan, University of Brest, France, or on a Flash EA1112 elemental analyser coupled to a Delta V Advantage isotope ratio mass spectrometer (Thermo Fisher scientific) at the LIENSs stable isotope facility, University of La Rochelle, France. Calibrations were based on reference materials (IAEA-600, IAEA-CH-6 and IAEA-N-2). Results were reported in the δ unit notation and expressed as parts per thousand (h) relative to the international standards (atmospheric N2 for nitrogen and Vienna-Pee Dee Belemnite for carbon). Analytical precision based on replicate measurements of acetanilide was < 0.15 h for both δ 13 C and δ 15 N values at both laboratories. Interlaboratory replicate measurements on reference materials (n=9), acetanilide (n=3), urea (n=3), casein (n=3), and lobster samples (n=3) found no significant differences in results between laboratories, with all differences below the analytical precision (<0.15 h). C:N ratios were determined from % element weights.

Fatty acid analysis

Subsamples of dry powders (n=141) were treated as described in [START_REF] Sardenne | Effects of extraction method and storage of dry tissue on marine lipids and fatty acids[END_REF].

Lipids were extracted from 50-100 mg of dry tissue with 6 mL of CHCl3:MeOH (2:1, v/v) mixture directly added into glass vials. Extracts were flushed with nitrogen, vortexed, sonicated for 15-20 min, and stored for 24 h at -20 °C. Tricosanoic acid (23:0) was added as an internal standard to 2 mL of lipid extract. Total lipids were transesterified with 800 μL of H2SO4 (3.8% in MeOH) at 100 °C for 10 min then washed with hexane-saturated distilled water. Fatty acid methyl esters (FAMEs) and dimethyl acetals (DMA; from vinyl ether lipids) were separated and quantified on a Varian CP8400 gas chromatography equipped with a Zebron ZB-WAX column (30 m length, 0.25 mm internal diameter, 0.25 μm film thickness; Phenomenex) and a flame ionisation detector at the Lipidocean facility, University of Brest, France. Samples were injected in splitless mode at 280 °C and carried by hydrogen gas. The oven temperature was raised from 60 °C to 150 °C at 50 °C min -1 , to 170 °C at 3.5 °C min -1 , to 185 °C at 1.5 °C min -1 , to 225 °C at 2.4 °C min -1 and then to 250 °C at 5.5 °C min -1 . FAMEs were identified by comparing sample retention times to those of four external standard mixtures (37-components, BAME, PUFA no. 1, and PUFA no. 3 FAME mix; Supelco) using Galaxie 1.9.3.2 software (Varian). The mean analytical variability of the method was 8.1%, based on Supelco 37-component FAME mix routinely checked. Individual FA results were expressed as a percentage of the total identified FAs, for 44 FAs accounting for > 1.5% in at least two samples. Total fat content was calculated as the sum of individual fatty acids (based on 23:0 concentration) and expressed in μg mg -1 of dry tissue.

2.4. Data analysis 2.4.1. Fractionation factors (FFs) and trophic level calculation Calculation of trophic levels and Bayesian mixing model include FFs (also called 'trophic enrichment factors' or 'trophic discriminant factors') to take into account for isotopic enrichment in heavier isotopes ( 13 C and 15 N) between successive trophic levels, i.e. between consumers and their prey. Although FFs are commonly set to 0.5h and 3.4h for δ 13 C and δ 15 N in marine systems, respectively, FFs are influenced by biological factors (e.g., species, tissue, age, food nutrients) [START_REF] Mccutchan | Variation in trophic shift for stable isotope ratios of carbon, nitrogen, and sulfur[END_REF][START_REF] Vander Zanden | Variation in δ 15 N and δ 13 C trophic fractionation: Implications for aquatic food web studies[END_REF] and a large range of FF values have been used to compute the trophic level of spiny lobsters based on δ 15 N values: e.g. 3.5h for Jasus edwardsii [START_REF] Guest | Trophic effects of fishing southern rock lobster Jasus edwardsii shown by combined fatty acid and stable isotope analyses[END_REF], 2.6h for P. cygnus (Waddington et al., 2008). To better take into account FFs' variability in the Seychelles marine food web, FFs were estimated for each potential prey of spiny lobster with simple linear regressions, as FFs tend to decrease with increasing isotopic ratios of prey [START_REF] Caut | Variation in discrimination factors (Δ 15 N and Δ 13 C): the effect of diet isotopic values and applications for diet reconstruction[END_REF]. FFs' estimation was thus based on the isotopic values of potential prey, and a diet experiment on P. cygnus [START_REF] Waddington | Diet quality and muscle tissue location influence consumerdiet discrimination in captive-reared rock lobsters (Panulirus cygnus)[END_REF]) (Fig. S1b). The obtained FFs were used to compute the trophic levels of spiny lobsters (see below) and were included in the Bayesian mixing model (see Section 2.4.2).

Trophic levels were calculated relative to the δ 15 N values of macroalgae (green, brown and red algae, average δ 15 N = 6.5 ± 0.5 h, n=18), according to the following equation:

𝑇𝑟𝑜𝑝ℎ𝑖𝑐 𝑙𝑒𝑣𝑒𝑙 = 𝛿 15 𝑁 𝑠𝑝𝑖𝑛𝑦 𝑙𝑜𝑏𝑠𝑡𝑒𝑟 -𝛿 15 𝑁 𝑚𝑎𝑐𝑟𝑜𝑎𝑙𝑔𝑎𝑒 𝐹𝐹 𝑚𝑎𝑐𝑟𝑜𝑎𝑙𝑔𝑎𝑒 + 1
with the average FF between lobster tail muscle and macroalgae previously computed from Fig.

S1b: FF15N = 2.9± 0.4 ‰). This calculation assumes the isotopic baseline to be the same for the three spiny lobster species.

2.3.2. Diet estimation using mixing modelling method A Bayesian mixing model was run to estimate the proportional contribution of sampled prey to the diet of spiny lobsters, based on SI values and FA profiles summarized through PCA. Because the discriminatory power of mixing models starts to decline above six or seven potential prey items [START_REF] Phillips | Best practices for use of stable isotope mixing models in food-web studies[END_REF], the prey species (n=13) were grouped as much as possible while respecting the phylogeny (n=9 classes) (Table 2). The model was run with two different groupings to compare (i) the diet among the three spiny lobster species, and (ii) the diet of spiny lobsters between reef types (granite vs carbonate reefs, regardless of the species because of low sample size and unbalanced design). The Bayesian mixing model included three types of parameters that were set as follow:

(i) Fractionation factors. For SI, FFs were obtained as described in Section 2.4.1. For FA profiles, FFs were set to 0 as the bioconversion of FA was assumed to be negligible in the hepatopancreas of spiny lobster due to its high proportion of dietary lipids [START_REF] Shu-Chien | Effect of dietary lipid source on expression of lipid metabolism genes and tissue lipid profile in juvenile spiny lobster Sagmariasus verreauxi[END_REF].

(ii) Concentration dependencies. Disparities among prey in the concentrations of trophic tracers induce different contributions of each prey to consumers' values [START_REF] Phillips | Incorporating concentration dependence in stable isotope mixing models[END_REF]. Percent of carbon and nitrogen elements and total fat content were thus included for each potential prey. Data analyses were performed using R software 3.5.0 (R Development Core Team, 2018), the JAGS program, and the ''simmr'' package (http://mcmc-jags.sourceforge.net) [START_REF] Parnell | Bayesian stable isotope mixing models[END_REF]. All results are reported as means ± SD.

Results

Variability in trophic tracers among potential lobster prey

Among the potential prey of spiny lobsters, the phylogenetic species class (i.e. 'Class') explained most of the variability in trophic tracers, while reef type explained no difference. δ 13 C values varied between species' class (F=39.6, d=8, p<0.001; 83.4% of explained variability), with values ranging from -23.8 ± 0.7 h for green algae of the Ulvophyceae class (Chlorodesmis spp.) to -10.8 ± 1.6 h for brown algae of the Phaeophyceae class (Turbinaria decurrens) (Fig. 2a). Reef type explained 0.5% of the variability in δ 13 C values (F=2.0, d=1, p=0.16). Similarly, δ 15 N values varied between species' class (F=59.2, d=8, p<0.001; 88.2% of explained variability), with values ranging from 6.0 ± 0.4 h for Turbinaria decurrens to 12.4 ± 0.1 h for the Echinoidea class (sea urchin Echinothrix diadema) (Fig. 2a). The reef type explained only 0.1% of the variability in δ 15 N values (F=0.3, d=1, p=0.59).

Regarding FA profiles, species' class and reef type explained 85.1% (χ2 = 645.9, d=8, p < 0.001) and 0.3% (χ2 = 1.9, d=1, p = 0.17) of the variability, respectively. PCA projection captured 39.5% the variability in FA profiles with 16:0, 18:0 DMA, 18:4n-1, 20:1n-11, 20:5n-3, and 22:5n-6 the six main discriminant FAs (Fig. 2b). PCA highlighted differences in FA profiles according to the phylogenetic groups: (i) the algae group (including Chlorophyta, Phaeophycea, and Rhodophyta phylums; PCA' bottom left panel) with the highest levels in 14:0, 16:0, 16:1n-7, and (F=2.2, df=2, p = 0.13; F=0.3, df=1, p = 0.58; and F=2.3, df=2, p = 0.11, respectively). However, δ 13 C values increased with carapace length (slope = 0.02; R² = 0.25; F=12.6, df=1, p < 0.001; 16.0% of explained variability), from -14.3 ± 0.6 h below 80 mm CL to -13.0 ± 0.3 h above 120 mm CL, especially for P. penicillatus (slope = 0.02, R² = 0.20; F=8.6, d=1, p < 0.01; Fig. 4).

Regarding trophic level, studied individuals ranged from 2.3 to 3.1, and only an interspecific difference was detected (F=32.1, d=2, p < 0.001; 51.9% of explained variability),i having a significantly lower trophic level than the two other species (mean trophic level of 2.4 vs. 2.8; Fig. 3b). No influence of carapace length, reef type, or interaction between species and carapace length on the trophic level was detected (F=0.1, d=1, p = 0.94; F=2.0, d=1, p = 0.17; and F=0.3, d=2, p = 0.75, respectively). Regarding FA profiles, only the carapace length was found to be an influencing factor (χ2 = 17.8, d=1, p < 0.001; 7.2% of explained variability). No influence of species, reef type, or interaction between species and carapace length on FA profiles was detected (χ2 = 4.1, d=2, p = 0.13; χ2 = 0.4, d=1, p = 0.55; and χ2 = 3.3, d=2, p = 0.19, respectively).

Model outputs

Inter-specific comparisons. Similar trends were observed for the three spiny lobster species, with Malacostraca identified as their main prey (average proportion: 47%-63% of the diet), followed by Echinoidea (16%-32% of the diet), Algae (8%-18% of the diet), Mollusc (5%-17% of the diet), and Anthozoa (1%-2% of the diet) (Fig. 5a). Due to large distributions in model outputs, dietary proportions of most potential prey were similar among spiny lobsters (at α = 0.05).

For instance, dietary proportions for Malacostraca were as follows: P. penicillatus 63 ± 14%, P. longipes 48 ± 10%, and P. versicolor 47 ± 12%, with 81.7% and 90.7% probability that P.

penicillatus consumed more Malacostraca than P. longipes and P. versicolor, respectively.

Noticeable differences mainly concerned P. versicolor: (i) 94.8% and 88.3% probability that P.

versicolor consumed less Echinoidea than P. longipes and P. penicillatus, respectively; (ii) 78.6% and 83.6% probability that P. versicolor consumed more Mollusc than P. longipes and P.

penicillatus, respectively, and (iii) 50.2% and 66.2% probability that P. versicolor consumed more Algae than P. longipes and P. penicillatus, respectively (Fig. 5a).

Reef habitats' comparison. Spiny lobsters caught from both carbonate reefs (n=7) and granite reefs (n=57) were found to mainly feed on Malacostraca (half of their diets; Fig. 5b).

However, spiny lobsters from granite reefs had 83.9% and 92.4% probability of feeding more on

Echinoidea and less on Molluscs than their congeners from carbonate reefs, respectively (Fig. 5b).

Discussion

Diet estimations from the mixing model were based on trophic tracers, and thus correspond to a diet assimilation over several days to months, based on turnover rates of SI and FA muscle and hepatopancreas [START_REF] Antonio | Influence of diet on the metabolic turnover of stable isotope ratios and fatty acids in the omnivorous shrimp Palaemon peringueyi[END_REF][START_REF] Waddington | Diet quality and muscle tissue location influence consumerdiet discrimination in captive-reared rock lobsters (Panulirus cygnus)[END_REF]. According to the mixing model, diet profiles of the three spiny lobster species were similar, with Malacostraca:

Anomura (hermit crabs) found to be their main prey. A similar finding has been observed from a reef lagoon of Puerto Morelos, Mexico, where hermit crabs and brachyurans were the main prey in the gut contents of P. argus juveniles (Briones-Fourzán et al., 2003). Preponderance of crustaceans in the diet of spiny lobsters has been observed in several Panulirus spp., usually in similar proportions to those of molluscs, depending on environmental factors such as reef type, season etc.

(Table S1). In our study, molluscs (gastropods + bivalves) ranked third to fourth in main prey items of spiny lobsters, which was lower than expected based on data from the literature (first to third main prey; Table S1). Around Mahé Island, bivalves do not form dense settlements (Selin et al., 1992;[START_REF] Taylor | Coral Reef and Associated Invertebrate Communities (Mainly Molluscan) Around Mahe, Seychelles[END_REF] and represent less than 14% of molluscs, which are dominated by gastropods (Barnes et al., 2009). It may explain why few bivalves were found by divers at the sampled sites (only the oyster Pinctada margaritifera). Consequently, bivalves contributed to a low proportion of the Seychelles spiny lobster' diets (<10% of the diets whatever the species), while they can be an important food source of Panulirus spp. in other regions of the world (e.g. Gulf of Mexico, Baja California; Table S1) and of other spiny lobsters e.g. Jasus spp. [START_REF] Mayfield | Diet of the West Coast rock lobster Jasus lalandii : influence of lobster size, sex, capture depth, latitude and moult stage[END_REF]Pollock, 1979). Echninoidea (sea urchins) were found to be an important diet component for the three species (16-32% of the diet). Our results are in accordance with observations from other regions, although Seychelles spiny lobsters seemed to rely more on Echinoidea compared to the 4-28% and 4-8% gut volume for P. guttatus and P. argus from Yucatan peninsula, Mexico, respectively [START_REF] Briones-Fourzán | Coral reef degradation differentially alters feeding ecology of co-occurring congeneric spiny lobsters[END_REF]. This relatively high consumption of sea urchin could be related to the high abundance of sea urchins in shallow areas around Mahé Island [START_REF] Pittman | Marine invertebrate communities of Baie Ternay National Marine Park and Baie Beau Vallon, Mahé, Seychelles. R Development Core Team[END_REF] and has already been suggested based on the high 20:4n-6 levels of both sea urchins and spiny lobsters in Seychelles (Sabino et al., 2020;[START_REF] Sardenne | Trophic structures in tropical marine ecosystems: a comparative investigation using three different ecological tracers[END_REF]. However, spiny lobsters from granite reefs seemed to feed more on Echinoidea and less on molluscs than individuals from carbonate reefs, while sea urchins are mainly found in patchy reefs of Mahé Island, i.e. reefs with a sand, rock or rubble base (Jan-Claas, 2020). This suggests spiny lobsters might move between reef types, despite narrow home ranges being suspected for these species and similar to results from related species (650 m² per week based on 50% utilisation distribution of P. interruptus from California, USA;

Withy-Allen and Hovel, 2013).

Interspecific differences in diet proportions were small, but P. versicolor consumed slightly more algae than the two other species, which is consistent with its lower trophic level (2.4 vs. 2.8).

This dietary difference could help to reduce interspecific competition, as even a small difference in trophic level favoured co-occurrence between Panulirus species on Caribbean reefs [START_REF] Segura-García | Dietary partitioning between sympatric species of spiny lobster in a coral reef system[END_REF]. Overall, the trophic levels of Seychelles spiny lobsters (ranging from 2.3 to 3.1, n=64) were similar to those derived from δ 15 N values for P. cygnus from Western Australia (1.9-2.2) (Waddington et al., 2008) and Jasus edwardsii from New Zealand (2.5-4.0) [START_REF] Jack | Individual variability in trophic position and diet of a marine omnivore is linked to kelp bed habitat[END_REF]. Such large ranges highlight the lobster's diet plasticity and their opportunistic feeding behaviour (Table S1). Regardless of the species, no ontogenetic change in trophic level was observed (i.e. no correlation between δ 15 N derived trophic level and carapace length). In contrast, δ 13 C values increased with spiny lobster size, indicating that large individuals may feed closer to the coast or more on prey with a strong benthic affinity (e.g. prey feeding on benthic organic matter with high δ 13 C values) than their smaller/younger congeners. While limited information is available on lobster size distribution with bathymetry in shallow areas, ontogenetic changes in dietary preferences are known for other spiny lobster species: e.g. it was observed for Jasus paulensis, J. lalandii and P. cygnus [START_REF] Blamey | Diet of the spiny lobster Jasus paulensis from the Tristan da Cunha archipelago: Comparisons between islands, depths and lobster sizes[END_REF][START_REF] Dumas | Strong preference for decapod prey by the western rock lobster Panulirus cygnus[END_REF][START_REF] Haley | Dietary change of the rock lobster Jasus lalandii after an "invasive" geographic shift: Effects of size, density and food availability[END_REF][START_REF] Mayfield | Diet of the West Coast rock lobster Jasus lalandii : influence of lobster size, sex, capture depth, latitude and moult stage[END_REF]. However, size related changes either in spatial distribution or in diet remain to be detangled.

The size ranges of the three species were too narrow in the present study to define size classes to be included in the mixing model, but it would be interesting to explore this aspect further by including juveniles in the dataset.

Study caveats

Caveats regarding mixing model outputs should be highlighted, since uncertainties in the estimated dietary proportions of algae, Malacostra, and Echinoidea remain strong despite a posteriori aggregations of specific prey classes, which could thus bias the results. These uncertainties have two main origins: (i) a low discrimination of some potential prey based on their SI compositions and FA profiles. Using additional trophic tracers such as sulphur or mercury stable isotopes might improve this discrimination (e.g. [START_REF] Higgs | Caribbean Spiny Lobster Fishery Is Underpinned by Trophic Subsidies from Chemosynthetic Primary Production[END_REF][START_REF] Jack | Prey base shifts in red rock lobster Jasus edwardsii in response to habitat conversion in Fiordland marine reserves:: implications for effective spatial management[END_REF]; and (ii) the influence of FF parameters, which were estimated from linear regressions based on one diet experiment only. Additional diet experiments would be required to accurately estimate FFs (for both SI and FAs) and thus increase robustness of the model setup as previously highlighted [START_REF] Brett | Using multiple tracers and directly accounting for trophic modification improves dietary mixing-model performance[END_REF][START_REF] Phillips | Best practices for use of stable isotope mixing models in food-web studies[END_REF].

Apart from the model parameters, more samples of spiny lobsters from carbonate reefs and of juveniles/small adults is needed. Due to loss of carbonate reefs in the region, samples of lobsters from this habitat were small and could prevent detection of any differences between reef types and ontogenetic stages. Some potential prey were also probably missing, especially bivalve and decapod species, which have reduced contributions to lobster's diets.

Future directions

Several integrated trophic indicators are particularly suitable for the implementation of an ecosystem-based fisheries management, including the monitoring of trophic levels of catches and the mean number of trophic links per species [START_REF] Tam | Towards ecosystem-based management: identifying operational food-web indicators for marine ecosystems[END_REF]. The present study will provide a baseline for trophic level of spiny lobsters from the Seychelles and might initiate monitoring. For instance, the seasonal monitoring of level would help to detect diet shifts in the past few months (e.g. following coral bleaching or after the fishing season). For the sampling to be nonlethal, muscle tissue can be sampled from the last segment of the fourth walking leg (ca. 0.5 g of wet tissue is required), which will start to grow again at the next molt, as after natural autotomy.

Monitoring may also allow trophic links per spiny lobsters species to be identified, especially regarding (i) potential shift in the diets of spiny lobsters before and after coral bleaching on carbonate reefs, as the dietary competition among spiny lobsters has increased after the coral bleaching event of 2016 (Sabino et al., 2020); and (ii) the influence of predator abundance, by comparing diets in heavily fished and protected regions throughout the ontogeny of spiny lobsters.

Spiny lobsters have several predators especially at juvenile stage, such as octopus, small sharks, groupers, and snappers, which influence their sheltering and foraging behaviours [START_REF] Ellis | Species interactions through ontogeny: Effects of size-selective predation by red grouper on Caribbean spiny lobster in solution holes of Florida Bay[END_REF][START_REF] Smith | Predation on early juvenile spiny lobsters Panulirus argus (Latreille): influence of size and shelter[END_REF].

Conclusion

Spiny lobsters from the Seychelles are generalist feeders that appear to mainly feed on crustaceans, sea urchins and algae. Dietary differences related to species, reef type and ontogeny were weak, but P. versicolor had a lower trophic level than the two other species. Lobsters from carbonate reefs fed more on Molluscs than those from granite reefs, and large lobsters probably fed closer to the coast or more on detritus feeders than the small ones. These results should however be refined by including more individuals (including size ranges) and prey types, trophic tracers, and parameters from diet experiments.

(

  iii) Priors. Model priors on diet source proportions should be included when potential sources have similar values of trophic tracers (Franco-Trecu et al., 2013). Priors on the nine prey class proportions were thus included, based on the literature on stomach contents for Panulirus spp. (

  18:3n-6; (ii) the Arthropod phylum (Malacostra class; PCA' top left panel), with the highest levels in 18:1n-7, 20:4n-6 and 20:5n-3; (iii) the Mollusc phylum (Bivalvia and Gastropoda classes; PCA' right panels), with the highest levels in three nonmethylated interrupted (NMI) FAs (22:2 NMI (7,13), 22:2 NMI (7,15), and 22:3 NMI), 16:1n-7 DMA, 22:1n-11, 22:5n-3, and 22:4n-6; (iv) the Echinoderm phylum (Echinoidea class; PCA' top right panel), with the highest levels in 18:0; and (v) the Cnidaria phylum (Anthozoa class), rich in 18:4n-1 (Fig. 2b;

  the mixing model did converge to 1. However, a posteriori groupings were required due to negative correlations between Malacostraca and Phaeophycea (Pearson's correlation coefficient r ranged from -0.57 to -0.92 according to lobster species and reef type), and to a lesser extent between Echinoidea and Malacostraca or Gastropoda. Despite the grouping of Chlorophyta, Ochrophyta, and Rhodophyta (under the 'Algae' group) and of Gastropoda and Bivalvia (under the 'Mollusc' group), a strong negative correlation remained for some groups (r of -0.87 between Malacostraca and Alga and of -0.69 between Echinoidea and Mollusc; Fig.S2), indicating the low discrimination of these groups in lobster diets and the possible overestimation of Malacostraca and Echinoidea proportions over Algae and Mollusc proportions in the spiny lobsters' diet.

Figure 1 .

 1 Figure 1. Sampling sites of spiny lobsters Panulirus longipes, P. penicilatus and P. versicolor

Figure 2 .

 2 Figure 2. Trophic tracers of spiny lobsters and their potential prey presented with (A) a biplot of nitrogen (δ 15 N) carbon (δ 13 N) isotopic values,

Figure 3 .

 3 Figure

Figure 4 .

 4 Figure 4. Outputs of the Bayesian mixing model estimating the prey proportions in the diet of spiny lobsters (A) among species, and (B) between habitat

  TableS1). To avoid this and to obtain a comprehensive overview of the lobsters' diet, proportions of prey classes were grouped a posteriori[START_REF] Phillips | Best practices for use of stable isotope mixing models in food-web studies[END_REF]) on a phylogenetic basis (i.e. Mollusc for Bivalvia and Gastropoda classes; Algae for Chlorophyta, Ochrophyta and Rhodophyta phylum).2.4.3. Variability in trophic tracers and diet estimationsFactors explaining variability in trophic tracers were assessed using ANOVA for univariate data (trophic level, δ 13 C and δ 15 N values), and PERMANOVA for multivariate FA profiles (nontransformed data). SIMPER analysis was used to detect FAs that contribute the most to differences among groups. For potential prey, the influence of the two fixed factors 'Class' (nine

	Samples for the posterior distributions were drawn by Markov Chain Monte Carlo (MCMC)

simulation methods (10,000 iterations, 4 parallel MCMC chains, and 1,000 burn-in rate). Model convergence diagnostics are presented in Fig.

S2

. Strong negative correlations between prey proportions may impair model estimations. phylogenetic classes) and 'Reef type' (granite vs carbonate reefs) was assessed. For spiny lobsters, the influence of 'Species' (P. longipes, P. penicillatus, and P. versicolor), 'Reef type' (granite vs carbonate reefs), 'Carapace length' and the interaction 'Species*Carapace length' were considered as fixed factors. Homogeneity of residuals was tested with Levene's test. Since homogeneity was found for all groups, parametric post hoc Tukey HSD tests were used to refine differences among groups. Diet estimations, i.e. outputs from the Bayesian mixing model, were compared among spiny lobster species and between reef types in the Bayesian framework (no test value, only probabilities were reported).

Table 2

 2

	).

Table 1 .

 1 Details on the sampling design of spiny losbters from reefs of Mahé Island, Seychelles.

	Species	Habitat type	N	Carapace length (mm)
	Panulirus longipes	Carbonate reefs	3	77.6 ± 5.1
	Panulirus longipes	Granitic reefs	23	78.7 ± 10.6
	Panulirus penicillatus	Carbonate reefs	4	119.1 ± 21.7
	Panulirus penicillatus	Granitic reefs	27	98.5 ± 16.0
	Panulirus versicolor	Granitic reefs	7	84.3 ± 15.5

Acknowledgments

This work is a contribution to the project SEYCCAT-BGF1/L/N6 "Science based restoration of commercially important spiny lobster habitats to help develop a sustainable fishery" led by the Marine Conservation Society Seychelles (MCSS); and to the Seychelles Participatory Lobster Monitoring Programme (PLMP) led by the Seychelles Fishing Authority (SFA). The project benefited from the financial support of the Seychelles' Conservation and Climate Adaptation Trust (SEYCCAT), the Seychelles Governement, the European Fisheries Partnership Agreement (EU-FPA), and the Institute for Research and Development (IRD). The authors would like to thank the lobster fishermen and the R.V. L'Amitie crew who assisted with sampling. The authors are also grateful to G. Guillou from the Plateforme Analyses Isotopiques of LIENSs, (University of La Rochelle, France) for his support in the lobster stable isotope analysis. Thanks are due to the CPER (Contrat de Projet Etat-Région) and the FEDER (Fonds Européen de Développement Régional) for funding the IRMS of LIENSs laboratory. The IUF (Institut Universitaire de France) is also acknowledged for its support to PB as a senior member.

Authors' contribution

Conceptualization: NB, LBa, LBl, DR; Funding acquisition: DR, NB, RG; Sampling: NB, LBa, LBl, KC, RM; Laboratory analyses: FS, SB, AB, JMR, FLL, FLG; Statistical analysis: FS; Writing


Julien Drevon

Florentin Millour

Pierre Cruzalebes

Peter Scicluna

Claudia Paladini

RHAPSODY: Reconstructing Hankel rAdial Profiles in centro-Symmetric Objects with Discrete rings for astrophYsics

Keywords: SPIE Proceedings, image reconstruction, intensity profile, astronomical objects

RHAPSODY is an intensity profile reconstruction code built to handle 2D centro-symmetric structures using interferometric data. RHAPSODY has been built to provide the community with a code requiring less parameters than classical image reconstruction code. This has as consequence leading to less solution non-uniqueness problems, a better convergence, and a better dynamic range assuming a centro-symmetric source at all wavelengths. RHAPSODY main steps are organized as following: In the first place, at each wavelengths, the code build a unique 1D structure made of concentric discrete uniform rings. Then, it applies the Hankel transform to reconstruct the equivalent visibility profile. Next, a change of coordinate is used on the 1D visibility profile to simulate the inclination and the rotation of the structure in the 2D Fourier plane. After a fitting process on the interferometric observations based on a χ 2 and Bayesian method, RHAPSODY apply an inverse Fourier transform and reconstruct the equivalent structure in the 2D image plane. According to preliminary tests made on RHAPSODY, the code is able to well reproduce 1D and 2D centrosymmetric structures up to a dynamic range of 0.5%.

INTRODUCTION

Image Reconstruction (IR) for long-baseline interferometry has been in the past a real challenge, that is still an active research topic today. An interferometer measures the amplitude and phase of the complex degree of coherence of the incoming light, meaning that in principle a straightforward inverse Fourier transform of interferometric data would enable one to get an image of the observed astrophysical object. However, this is not possible in practice because of several issues: an interferometer provides an uneven sampling of the data, leading to information voids ; it provides a non-linear way of making available the measures (amplitude and phase), leading to a non-convex problem ; finally an image contains a number of information (pixel intensities) that outperform by a large fraction the limited number of measurements, leading to a non-unique solution. Solving the IR problem in interferometry is therefore solving a non-convex inverse problem in a sparse case. Thanks to innovative ideas like the use of a Bayesian approach combined with maximum entropy 1 or various other regularizations, 2 plus an efficient weighting and grid-search, [START_REF] Hofmann | An image reconstruction method (IRBis) for optical infrared interferometry[END_REF] self-calibration, [START_REF] Millour | Hot Stars and Interferometry[END_REF] and more recently the emergence of machine learning, [START_REF] Claes | Neural network based image reconstruction with astrophysical priors[END_REF] IR has become possible with long-baseline interferometry.

Yet, IR is still of a difficult access for the general user, because mainly of the non-unicity of the solution, and the hand-tweaking of parameters to solve the IR problem (albeit it has been made available, e.g. with oi-imaging). One way to tackle these issues is to reduce by a large fraction the number of input parameters to the inverse problem, either by considering the imaging problem intrinsically sparse, [START_REF] Sanchez-Bermudez | The 2016 interferometric imaging beauty contest[END_REF] or by radically changing the IR paradigm.

RHAPSODY is a proposition for the second assertion: instead of reconstructing an image, this tool proposes to reconstruct radial profiles.

Further author information: Julien Drevon: E-mail: julien.drevon@oca.eu, Telephone: +33 (0)6 11 93 42 86 1 This makes it a very complementary tool to already existing image reconstruction software. The aims of this code is to provide a fast and robust geometrical description of an astrophysical environment with a centrosymmetric model. We adopt a 1D Hankel transform combined with a Bayesian approach to reconstruct the radial intensity profile of the source. The Hankel transform has the advantage of requiring less parameters than classical image reconstruction techniques, reducing the non-unicity of solutions, improving convergence, and providing a better dynamic range.

In Section 2 we present step by step the RHAPSODY tool and its features to get an insight in the core of the code and the process involved in the "image" reconstruction. Next, in Section 3 we discuss the robustness of the RHAPSODY outputs thanks to simulated interferometric data and we estimate its intrinsic dynamic range. Then, finally we show in Section 4 some examples of RHAPSODY applied to real or simulated interferometric data.

THE RHAPSODY TOOL AND FEATURES

RHAPSODY * can handle interferometric data that are stored in the standard OIFITS format. [START_REF] Duvert | OIFITS 2: the 2nd version of the data exchange standard for optical interferometry[END_REF] In this section we will go through all the important steps of the software.

Setting the initial parameters

The important initial parameters are only a few. We can separate them in 3 different categories: i) facilitate the work of reading OIFITS files, ii) init the model structure, iii) define the fitting method.

To facilitate the work of reading OIFITS files, RHAPSODY needs as input the number of telescope involved in the observations. Furthermore, in order to avoid the presence of unexpected correlated flux among the OI VIS and OI VIS2 table, the user needs to precise if the FLUX table is expected in the data. If not, we assume that in absence of the FLUX table, the OI VIS and OI VIS2 tables contain visibilities and squared visibilities, respectively. In the case where several bands are present in the OIFITS file, the user need to provide the names and the wavelengths intervals in which he wants to do the image reconstruction.

To initialize the model structure, RHAPSODY needs the initial intensity profile as well as the number of rings and their respective angular diameter (see Section 2.2 for details). Up to now, the user can only chose among the Gaussian, Uniform Disk, or power-law profiles as initial guess.

Finally, to define the fitting method, RHAPSODY needs the following inputs: the desired regularization (between the Total Variation and the Quadratic Smoothness, see Section 2.3 for details), the value of the hyperparameter to use (a list of hyper-parameter values can be provided), the fitting method (by default COBYLA), as well as the tolerance and the maximum of iterations to converge to a solution (the computation time is very sensitive to the tolerance).

Rings and Hankel Profile

Instead of forcing the routine to describe a given component (for example the dust envelope) with a particular geometrical model (for example a Gaussian disk), we use a large number of concentric uniform narrow rings (of width given as an input) to reproduce the radial profile of the stellar environment. With this method, thick rings are reproduced by contiguous narrow rings. This Hankel profile is built thanks to a set of N step concentric circular rings with diameters θ k = kθ 0 , where θ 0 = 0.5λ/B being selected as a fraction of the angular resolution of the interferometer. Within that assumption, we write the normalized monochromatic radial profile as:

I Hankel [h k (λ)] = 4 π Nstep-1 k=0 h k (λ) Π(θ k+1 ) -Π(θ k ) θ 2 k+1 -θ 2 k (1)
where Π(θ) is the rectangular function of width θ. The height in intensity h k of each step fulfills the condition 0 ≤ h k ≤ h 0 whatever the wavelength is. This is to avoid solutions having a surrounding environment brighter than a central star.

Thus, the final expression of the global visibility is:

V Hankel [q, a k (λ)] = Nstep-1 k=0 a k (λ) θ 2 k+1 V UD (q; θ k+1 ) -θ 2 k V UD (q; θ k ) θ 2 k+1 -θ 2 k , (2)
where the visibility coefficients a k , that are varied by our reconstruction method, are:

a k (λ) = h k (λ) Nstep-1 k=0 h k (λ) (3)
One can note here the little trick of using uniform disks to represent rings with a non-zero thickness. The visibilities of each uniform disk can be pre-computed initially to save computing time afterwards. A consequence of this trick is the step-like appearance of our radial profiles. A smoother version is in preparation and its description is beyond the scope of this paper.

Adjusting V 2

Hankel to the measured squared visibilities V 2 obs (with their uncertainties σ obs ) is then a matter of minimizing the distance between both through a χ 2 minimization, varying the a k coefficients.

χ 2 (λ) = q V 2 Hankel [q, a k (λ)] -V 2 obs [q] 2 σ 2 obs (4)
However, trying to adjust the visibility deduced from a lot of regular intensity steps using a bare final χ 2 value usually gives poor results, because the number of parameters a k to fit becomes large compared to the number of available data, in a similar way as in an image reconstruction problem. This is where the Bayesian method comes into play by adding to the simple χ 2 an a priori information.

The Bayesian Approach

This additional information to the problem-solving is called a regularization term f prior , calculated on the variable parameters themselves. Such a regularization term must be designed so as to maximize the probability that the reconstructed profile resembles the "true" one. For more details on how to deal with heavy inverse problems in a Bayesian framework, we refer the reader to the excellent introduction to image reconstruction for optical interferometry. [START_REF] Thiebaut | Image reconstruction in optical interferometry[END_REF] Adjusting the Hankel profile to the visibility data is then a matter of adding the χ 2 and the regularization term f prior multiplied with a so-called hyper-parameter µ, that is arbitrarily set:

f tot (λ) = χ 2 (λ) + µ f prior (λ) (5)
In our case, we have tested two different regularization terms for f prior : (1) a smoothness term f prior = f smooth , or (2) a total variation term f prior = f var . On the one hand, the smoothness f smooth , described by the intensity difference between adjacent rings, is given by:

f smooth (λ) = Nstep-1 k=0 h k+1 (λ) -h k (λ) 2 . ( 6
)
f smooth will tend to minimize the variations between consecutive radii, resulting into a smooth radial profile. On the other hand, the total variation regularization term f var defined as:

f var (λ) = Nstep-1 k=0 h k+1 (λ) -h k (λ) θ k+1 -θ k 2 , (7)
will tend to reconstruct uniform intensities with sharp edges. Solving the problem (i.e. finding the right radial profile) consists therefore in minimizing the term in Eq. 5, by varying the parameters a k using a gradient descent or any other optimization algorithm (in our case: the Cobyla minimization routine from the Lmfit package † library that can be found in Python). This description of the algorithm corresponds exactly to the version of RHAPSODY the has been been used in [START_REF] Drevon | Locating dust and molecules in the inner circumstellar environment of R Sculptoris with MATISSE[END_REF]

Rotation and Inclination

In order to make progress towards a more versatile image reconstruction tool, we introduced at the occasion of the 2022 interferometric imaging beauty contest the possibility to reconstruct flattened structures (like e.g. disks).

All the transformations are done directly in the Fourier plane with coordinates (u, v). The inclination i is simulated by a cos(i) compression factor along the chosen axis which becomes the major-axis in the Fourier plane (i.e., the minor axis in the image plane).

Fourier plane :

u ′ = u v ′ = v/ cos(i) Image plane : x ′ = x y ′ = y * cos(i) (8)
Then, the rotation is made through a change of coordinate with a position angle θ:

Fourier plane :

u " = u ′ sin(θ) + v ′ cos(θ) v " = -u ′ cos(θ) + v ′ sin(θ)
Image plane :

x " = x ′ cos(θ) + y ′ sin(θ)

y " = -x ′ sin(θ) + y ′ cos(θ) (9)

From 1D to 2D Fourier plane and "Image" Reconstruction

To reconstruct a flattened radial profile, we first need to pick the pre-computed visibilities at the right compressed and rotated coordinates. Since before the change of coordinate, we are dealing with angular centro-symmetric structure, all the information which is shrunk by the compression factor and rotated by an angle θ is located inside the 1D visibility curve. Hence, RHAPSODY makes an interpolation at spatial frequencies q θ = u " 2 + v " 2 on the 1D visibility curve to get the new visibility points of the inclined and rotated modeled structure. Once we have a complete 2D Fourier plane, the code apply an inverse Fourier transform and retrieve the equivalent structure in the image plane.

OUTPUT ROBUSTNESS

Here we test RHAPSODY against simple and more complex models to assess its intrinsic capabilities and limits.

Uniform Disk

We modeled interferometric data based on multiple uniform disks with different diameter using the ASPRO2 software of the JMMC ‡ . We simulated these uniform disks observations as if they were seen by the VLTI/MATISSE instrument. The aim is double: first to check that RHAPSODY is able to correctly reproduce simple model, and second to reach the limit of the code for a given interferometric data-set. In those simulations, we have selected a hyperparameter of µ = 1. The resolution of the MATISSE instrument in L-band is ranging between 3.6 and 275 mas. The results are shown in figures 1, 2, and 3.

First, we observe that at large scales in Fig. 1, RHAPSODY is able to reproduce the uniform disks up to a ≈ 100 mas radius, the 250 mas uniform disk is not well reproduced by RHAPSODY as we reach the interferometric field of view limit of the instrument.

Figure 2 highlights lower scales that the instrument is able to resolve, we can see that RHAPSODY does reproduce the boundaries of the uniform disk properly. We note the presence of artifacts at low flux, especially for the 10 and 20 mas radius uniform disks.

The last modeling attempt is at lower scales (Figure 3), we can see that even if the theoretical resolution of the instrument is at 3.6 mas, we are still able to reproduce structures up to at least half the resolution of the ‡ made freely available at http://jmmc.fr/ ~webmaster/jmmc-html/aspro_page.htm instrument. However we remark that the boundaries of the structures are well reproduced, but the shape of the structure is closer to a Gaussian shape than of an uniform disk.

Hence, we can see here the limits of RHAPSODY: At large scales (at or above the interferometric field of view λ/B min ), the code struggle to proper produce the exact shape of the object, and at small scales, RHAPSODY is still able to get the correct size of the structure up to about half the resolution of the instrument, in a very similar way as for image reconstruction. [START_REF] Thiebaut | Image reconstruction in optical interferometry[END_REF]

Gaussian Disk

We made a test with an Gaussian model with a 10 mas FWHM in the same way as for uniform disks (Figure 4). We can assert here that a Gaussian model is also generally well reproduced by RHAPSODY.

However, the dynamic range is limited to about 100, as can be seen at the end of the Gaussian wings: at intensities levels between 1% and 0.1% of the peak intensity, the faint reconstructed region starts to exhibit similar artifacts as in Figure 2, hence the flux starts to be poorly spatially constrained.

Elongated Gaussian Ellipse

In this subsection we wanted to test the new feature which bring the 2D aspect of RHAPSODY involving inclination and rotation. We took the same Gaussian disk as in section 3.2 but we added an inclination of 60 • and a position angle of 80 • and we used RHAPSODY to reconstruct the environment using the modeled interferometric data with ASPRO. The method involved in RHAPSODY allows us to have a face-on description of the object in addition to the reconstructed image produced at the end. As seen in Fig 5, we have successfully reconstructed the Elongated Gaussian Ellipse using RHAPSODY and a hyperparameter of µ = 10 3 . The output are the same as the input of our model: as much in terms of the Gaussian FWHM as well as in terms of inclination angle and rotation, with a precision better than 1%.

1D Dust Radiative Transfer environment

In this paragraph we test the robustness of the RHAPSODY against a realistic model used in our paper on R Scl. [START_REF] Drevon | Locating dust and molecules in the inner circumstellar environment of R Sculptoris with MATISSE[END_REF] We built a toy model mimicking the R Scl environment using DUSTY and a Gaussian layer at 15 mas with an intensity ratio with respect to the central source starting from 1/1000 and reaching up to 1/10 (Fig 8 in appendix). Using ASPRO2, we generate VLTI/MATISSE visibilities, add error noise to the data, and use the same (u,v)plane coverage for each of our different models and for both LM -and N -bands.

The initial guess intensity profiles used by RHAPSODY for all the simulations in both bands are similar to the one used in this paper (i.e. Gaussian profile with a F W HM = 7 mas). The value of the hyperparameter used is µ = 1. Figure . 8 shows the spectra of the intensity profile reconstruction produced by RHAPSODY. We can notice that the lower limit of the dynamic range in L-band is about 1/200 (0.5%). In N -band, the code does not reconstruct (or very marginally) the inner Gaussian structure added at 15 mas whatever its intensity. We suspect this non-detection in N -band is due to an unfortunate combination of instrumental angular resolution and dynamic range.

USING RHAPSODY ON REAL AND BENCHMARK DATA

4.1 Spherical envelope: R Sculptoris R Scl was observed for the first time simultaneously in both L-and N -bands with the VLTI/MATISSE instrument. We used RHAPSODY to analyze this dataset. [START_REF] Drevon | Locating dust and molecules in the inner circumstellar environment of R Sculptoris with MATISSE[END_REF] Thanks to this tool, we have discovered interesting structures in the vicinity of the star. As we can see in Fig 6, RHAPSODY has spotted two interesting structures in L-band: one between 5 and 7 mas, above the Rosseland radius (dashed blue line) and one between the two green lines between 12 and 18 mas. A layer beyond the estimated Rosseland radius of R Scl is seen as a faint environment between 5 and 7 mas (∼ 1.0-1.3 R ⋆ ), assuming that the border of the photosphere is at 5 mas from the center of the star.

Thanks to the spectral information provided by the MATISSE instrument we interpret this layer as a mixture of molecules (C 2 H 2 + HCN). In addition, we have spotted another layer at a further distance from the star which could be made of both: dust and molecules.

Inclined Disk: SPIE Beauty Contest 2022

As seen in Figure 7, we also took the opportunity of this SPIE 2022 conference to participate to the SPIE Beauty Contest 2022. We used a total of 56 rings with a 2 mas width, we have a total of 616 parameters for the 11 wavelengths involved in RHAPSODY. We added two parameters in our algorithm: flattening and position angle on the plane of the sky, for each wavelengths, leading to a total of 638 parameters for 1440 visibility points. Hence, the number of parameters is inferior to the number of measurements.

We applied it to the first set of the beauty contest "contest1". The RHAPSODY images show a bright central source surrounded by an elongated 8.5 mas radius (major axis) ring with a flattening ranging from 0.6 to 0.7. We can compute an inclination of that ring, that ranges from 45 to 50 • , corresponding to the inclination angle of a disk if that structure were indeed a thin disk. The disk position angle ranges from 41 to 55 • . The central source seems to be resolved or partially resolved.

Figure 1 :

 1 Figure 1: Left Panel: Intensity profile spectra of uniform disk modeled data of respectively a diameter of 100, 200, 500 and 700 mas (50, 100, 250 and 350 mas radius). Right panel: RHAPSODY intensity profile spectra outputs.

Figure 2 :

 2 Figure 2: Left Panel: Intensity profile spectra of uniform disk modeled data of respectively a diameter of 4, 10, 20 and 40 mas (2, 5, 10 and 20 mas radius). Right panel: RHAPSODY intensity profile spectra outputs.

Figure 3 :

 3 Figure 3: Left Panel: Intensity profile spectra of uniform disk modeled data of respectively a diameter of 0.2, 0.4, 1 and 2 mas (0.1, 0.2, 0.5 and 1 mas radius). Right panel: RHAPSODY intensity profile spectra outputs.

Figure 4 :

 4 Figure 4: Left Panel: Intensity profile spectra of modeled data of a Gaussian with a FWHM of 10 mas. Right panel: RHAPSODY intensity profile spectra of the same modeled data.

Figure 5 :

 5 Figure 5: Left Panel: RHAPSODY simulations results of the Elongated Gaussian Ellipse with a FWHM of 10 mas, face-on. Right panel: Same RHAPSODY results but as seen by the observer.

Figure 6 :

 6 Figure 6: Upper panels correspond to the observed MATISSE spectrum normalized by the median and then divided by a black body spectra at T=2700 K to underline the observed emission and absorption features. Lower panels are the spectro-radial maps in L-and N -band obtained by plotting the best intensity Hankel profile normalized at one for each observed wavelength. The faint red structures highlighted with the inclined dashedwhite lines are reconstruction artifacts). The dashed and solid blue horizontal lines show the position of the Rosseland radius and the DUSTY inner radius respectively. The solid green horizontal lines delimit the extension of a hot distinct molecular layer). The grey vertical bands cover the spectral ranges where the centro-symmetric Hankel profile is not able to properly reproduce the asymmetric shape of R Scl revealed by non-zero closure phases.

Figure 7 :

 7 Figure 7: RHAPSODY contribution to the SPIE 2022 beauty contest.

  (c) Additional gaussian layer at 15 mas with an intensity ratio with respect to the star of 1/200 (d) Additional gaussian layer at 15 mas with an intensity ratio with respect to the star of 1/100

Figure 8 :

 8 Figure 8: Top view continued...

* made freely available at https://github.com/jdrevon/RHAPSODY

† made freely available at https://github.com/lmfit/lmfit-py/

§ Link to ASPRO on www.jmmc.fr

ACKNOWLEDGMENTS

MATISSE has been built by a consortium composed of French (INSU-CNRS in Paris and OCA in Nice), German (MPIA, MPIfR and University of Kiel), Dutch (NOVA and University of Leiden), and Austrian (University of Vienna) institutes. It was defined, funded and built in close collaboration with ESO. The Conseil Départemental des Alpes-Maritimes in France, the Konkoly Observatory and Cologne University have also provided resources to manufacture the instrument. The authors made extensive use of the Jean-Marie Mariotti Center (JMMC)