Guillaume Marchessaux
email: guillaume.gmarchessaux@gmail.com

Mireille Harmelin-Vivien

Mélanie Ourgaud

Daniela Bănaru

Loïc Guilloux

Bruno Belloni

Benoit Lebreton

Gaël Guillou

Delphine Thibault

First overview on trophic relationships of the invasive ctenophore Mnemiopsis leidyi in a Mediterranean coastal lagoon (Berre Lagoon, France): benthic-pelagic coupling evidenced by carbon and nitrogen stable isotope composition

Keywords:

Introduction

Coastal lagoons play a major role in coastal biodiversity and human societies supporting high biological productivity, economic and socio-cultural activities, etc. [START_REF] Lemée | What meaning do individuals give to coastal risks? Contribution of the social representation theory[END_REF]. However, they are very sensitive to environmental changes such as nutrient over-enrichment, contamination, habitat destruction, and non-indigenous species introductions [START_REF] Audouit | Comparing social representation of water quality in coastal lagoons with normative use of ecological indicators[END_REF]. Non-indigenous species (NIS) introductions in new ecosystems constitute a major source of biological pollution [START_REF] Elliott | Biological pollutants and biological pollution--an increasing cause for concern[END_REF][START_REF] Boudouresque | An overview of species introduction and invasion processes in marine and coastal lagoon habitats[END_REF] and are often considered as invasive strong ecological and economic impacts, altering the diversity and functioning of ecosystems, and negatively impacting human activities such as fisheries, industry and tourism [START_REF] Streftaris | Globalisation in marine ecosystems: the story of non-indigenous marine species across European seas[END_REF][START_REF] Gallardo | Global ecological impacts of invasive species in aquatic ecosystems[END_REF]. Understanding invasive * Corresponding author.

NIS impacts on ecosystem functioning requires the determination of their position and trophic relationships within the food web [START_REF] Byrnes | Invasions and extinctions reshape coastal marine food webs[END_REF]. The ctenophore Mnemiopsis leidyi Agassiz 1865 is ranked among the 100 worst bio-invaders in the world by IUCN [START_REF] Lowe | The Invasive Species Specialist Group (ISSG) a specialist group of the Species Survival Commission (SSC) of the World Conservation Union (IUCN), p. 12[END_REF], with the source populations located in estuaries and coastal waters of the eastern seaboard of the Americas, from Chesapeake Bay on the east coast of the USA through the Gulf of Mexico to the coast of Argentina [START_REF] Gesamp | Opportunistic settlers and the problem of the ctenophore Mnemiopsis leidyi invasion in the Black Sea[END_REF]. This species has spread for the past 30 years through Eurasia following two main pathways [START_REF] Ghabooli | Multiple introductions and invasion pathways for the invasive ctenophore Mnemiopsis leidyi in Eurasia[END_REF][START_REF] Ghabooli | Invasion Pathway of the Ctenophore Mnemiopsis leidyi in the Mediterranean Sea[END_REF][START_REF] Bayha | Worldwide phylogeography of the invasive ctenophore Mnemiopsis leidyi (Ctenophora) based on nuclear and mitochondrial DNA data[END_REF]. The first one connected the Gulf of Mexico seeding population to the Black Sea where in 1982 M. leidyi was first reported [START_REF] Pereladov | Some observations for biota of Sudak Bay of the Black Sea[END_REF]. It then spread rapidly to the Caspian Sea, and to various waters of the Mediterranean basin (see Shiganova et al., 2019, for a review). The second invasion linked the east coast of the USA to the North Sea and the Baltic Sea [START_REF] Javidpour | First record of Mnemiopsis leidyi a. Agassiz 1865 in the Baltic Sea[END_REF][START_REF] Oliveira | The presence of the ctenophore Mnemiopsis leidyi in the Oslofjorden and considerations on the initial invasion pathways to the North and Baltic Seas[END_REF][START_REF] Van Ginderdeuren | Distribution of the invasive ctenophore Mnemiopsis leidyi in the Belgian part of the North Sea[END_REF][START_REF] Antajan | The invasive ctenophore Mnemiopsis leidyi A. Agassiz, 1865 along the English Channel and the North Sea French coasts: another introduction pathway in northern European waters?[END_REF]. Along the French Mediterranean coasts, M. leidyi was reported from early 2000 in several coastal lagoons, from Bages-Sigean in the west to Berre Lagoon in the east [START_REF] Delpy | Identifying the drivers of abundance and size of the invasive ctenophore Mnemiopsis leidyi in Northwestern Mediterranean lagoons[END_REF][START_REF] Marchessaux | Etude de L'Anthroposystème Emblématique de L'étang de Berre : Approches écosystémique Et Sociologique de L'Impact Du Cténaire Invasif Mnemiopsis Leidyi[END_REF] and in Corsica [START_REF] Shiganova | Patterns of invasive ctenophore Mnemiopsis leidyi distribution and variability in different recipient environments of the Eurasian seas: A review[END_REF].

The arrival and establishment of a large population of M. leidyi disturbed local environmental conditions both at the ecological (high predation pressure, competition with native planktivores) and economic level (i.e. clogging nets, catch reduction; [START_REF] Vinogradov | Consumption of zooplankton by the comb jelly Mnemiopsis leidyi and pelagic fishes in the Black Sea[END_REF]Purcell et al., 2001;[START_REF] Marchessaux | Etude de L'Anthroposystème Emblématique de L'étang de Berre : Approches écosystémique Et Sociologique de L'Impact Du Cténaire Invasif Mnemiopsis Leidyi[END_REF]. The diet of M. leidyi depends on different parameters such as ontogeny, size, food availability and season. Adults have a broad meso-zooplanktivorous diet, including fish eggs and larvae [START_REF] Purcell | Interactions of pelagic cnidarians and ctenophores with fish: a review[END_REF][START_REF] Purcell | Extension of methods for jellyfish and ctenophore trophic ecology to large-scale research[END_REF], while larvae feed on microplankton [START_REF] Rapoza | Life-stage dependent, in situ dietary patterns of the lobate ctenophore Mnemiopsis leidyi Agassiz 1865[END_REF][START_REF] Sullivan | Gut evacuation of larval Mnemiopsis leidyi A. Agassiz (Ctenophora, Lobata)[END_REF].

Gelatinous species in general (pelagic ctenophores and cnidarians) are often grouped into a single trophic category in ecosystem models [START_REF] Pauly | Jellyfish in ecosystems, online databases, and ecosystem models[END_REF][START_REF] Condon | Questioning the rise of gelatinous zooplankton in the world's oceans[END_REF] and are considered as direct competitors of planktivorous fish [START_REF] Sommer | Pelagic food web configurations at different levels of nutrient richness and their implications for the ratio fish production: primary production[END_REF][START_REF] Brodeur | Spatial overlap and trophic interactions between pelagic fish and large jellyfish in the northern California Current[END_REF]. However, they cover a wide range of species and diets, including predators of other gelatinous zooplankton [START_REF] Greve | Plankton time-space dynamics in German Bighta systems approach[END_REF]. It is therefore essential to assess trophic diversity at the species level [START_REF] Nagata | Food web characterization based on δ 15 N and δ 13 C reveals isotopic niche partitioning between fish and jellyfish in a relatively pristine ecosystem[END_REF]. Several techniques have been used to study the trophic ecology of gelatinous zooplankton, such as traditional gastrovascular cavity content analyses [START_REF] Larson | In situ feeding rates of the ctenophore Mnemiopsis mccradyi[END_REF][START_REF] Hansson | A method for in situ estimation of prey selectivity and predation rate in large plankton, exemplified with the jellyfish Aurelia aurita (L.)[END_REF] and grazing experiments [START_REF] Mazlum | Gastric emptying, clearance rate, feeding periodicity and food consumption of the Black Sea jelly fish, Mnemiopsis leidyi (Agassiz)[END_REF][START_REF] Rowshantabari | Effect of temperature on clearance rate, daily ration and digestion time of Mnemiopsis leidyi from the southern Caspian Sea[END_REF], which indicate which prey was lately ingested, and more recently biochemical tracers, such as stable isotopes (SI) and fatty acids, which provide an integrated view of the assimilated food sources [START_REF] Pitt | Stable isotope and fatty acid tracers in energy and nutrient studies of jellyfish: a review[END_REF][START_REF] Vansteenbrugge | Trophic ecology of Mnemiopsis leidyi in the southern North Sea: a biomarker approach[END_REF]. Fractionation from one trophic level to the next is assumed to be generally from 0.1 to 1 for δ 13 C and from 2 to 4 for δ 15 N, which allows the use of stable isotopes to determine the trophic position of organisms within the food web (e.g. [START_REF] Peterson | Stable isotopes in ecosystem studies[END_REF][START_REF] Layman | Applying stable isotopes to examine food-web structure: an overview of analytical tools[END_REF][START_REF] Glibert | Stable isotope tracers: Enriching our perspectives and questions on sources, fates, rates, and pathways of major elements in aquatic systems[END_REF]. They also provide a way to quantify the contributions of different food sources if they statistically differ [START_REF] Parnell | Source partitioning using stable isotopes: Coping with too much variation[END_REF][START_REF] Phillips | Best practices for use of stable isotope mixing models in food-web studies[END_REF]. SI composition of prey can identify common resources, which can provide a basis for hypotheses regarding competition or predation interactions among consumer organisms [START_REF] Hamer | The invasive ctenophore Mnemiopsis leidyi: a threat to fish recruitment in the North Sea?[END_REF][START_REF] Frost | Distribution and trophic links of gelatinous zooplankton on Dogger Bank, North Sea[END_REF][START_REF] Kellnreitner | Feeding interactions between the introduced ctenophore Mnemiopsis leidyi and juvenile herring Clupea harengus in the Wadden Sea[END_REF]. A recent modeling study of the growth of M. leidyi indicates that variations in food quantity and quality may induce rapid changes in its population dynamics [START_REF] Alekseenko | Main factors favoring Mnemiopsis leidyi individuals growth and population outbeaks: A modelling approach[END_REF][START_REF] Marchessaux | Refugia area for the ctenophore Mnemiopsis leidyi A. Agassiz 1865 in the Berre Lagoon (southeast France): The key to its persistence[END_REF], but only planktonic prey were included in the biogeochemical model used. However, some studies demonstrate that, in coastal environments, jellyfish and ctenophores are able to prey on benthic organisms when they emerge from the sediment at night [START_REF] Pitt | Predation by jellyfish on large and emergent zooplankton: Implications for benthic-pelagic coupling[END_REF] or benthic organic matter [START_REF] Cui | Potential dietary influence on the stable isotopes and fatty acid compositions of jellyfishes in the Yellow Sea[END_REF][START_REF] Javidpour | Temporal dietary shift in jellyfish revealed by stable isotope analysis[END_REF], playing therefore a role in the benthic-pelagic coupling.

The present study is a first attempt at elucidating the feeding relationships of Mnemiopsis leidyi with planktonic and benthic compartments in Berre Lagoon (Mediterranean Sea) through a stable isotope approach, which was never explored before in this coastal lagoon. Its aims were to (1) provide the C and N stable isotope composition of M. leidyi and its main food components, (2) evidence the variations of C and N isotope composition with individual size, (3) quantify the respective contributions of pelagic and benthic sources of prey in its diet, and (4) consider the eventual trophic relationships between M. leidyi and the jellyfish Aurelia sp.

Materials and methods

Study area

Berre Lagoon is the largest (surface area 155 km 2 , total volume of 0.98 x10 9 m 3 , average depth 7 m) semi-enclosed brackish lagoon located near Marseille, southeast France. Inflows of freshwater into the lagoon are all situated in the northern part with two natural rivers (Touloubre and Arc), and the wide industrial channel of the Durance supplying several hydroelectric powerplants, the nearest located at Saint-Chamas (Fig. 1). These freshwater inputs participate to the high variability of the environmental parameters of the Berre Lagoon ecosystem, where temperature and salinity in particular largely fluctuate over the year (Table 1). Berre Lagoon is connected to the Mediterranean Sea in its southwestern part through the Caronte channel. Low-amplitude (∼30 cm) semi-diurnal tides are observed. Berre Lagoon is composed of two sub-basins: the larger basin (∼9 m maximum depth) and the Vaine basin (∼5 m maximum depth), separated by a sand shoal (∼1.5 m depth). Extensive Zostera noltii meadows are observed particularly in the Vaine basin.

Sampling and sample preparation

Sampling was done at a station located in the Vaine basin (latitude: 43.446, longitude: 5.187; depth: 4 m) (Fig. 1), on 20 April 2016. Samples for benthic (sediment organic matter -SOM, macrophytes, invertebrates and juveniles of teleostean fish) and for pelagic (particulate organic matter -POM, zooplankton, the ctenophore Mnemiopsis leidyi, and the cnidarian Aurelia sp.) compartments were collected by divers or nets.

Benthic compartments: SOM was collected by divers by scraping the surface sediment with a sterile 1 L plastic bottle (∼500 g of matter). Four independent replicates were performed. In the laboratory, bottle contents were sieved on a 200 µm mesh sieve to remove large sediment particles and kept in polyethylene bags. Several individuals of the seagrass Zostera noltii and the green alga Chaetomorpha linum were collected by divers. In the laboratory, rhizomes, and damaged parts of the macrophytes were discarded and epiphytes gently removed by scraping the leaves. The clean leaves were rinsed with ultrapure water and placed in polyethylene bags.

Pelagic compartment: Four water samples (POM) were collected ∼1 m below the surface using a sampling bottle. Between 2 to 4 L of water were first sieved through a 80 µm mesh silk to remove detritus and zooplankton organisms, and filtered on pre-combusted (450 • C for 5 h) and pre-weighed GF/F filters (porosity: 0.7 µm). Then filters were stored in hermetically sealed boxes. Zooplankton was collected using a vertical 80 µm mesh net (1.2 m long, 30 cm opening area), the content of the codend was then sieved onto a column of three sieves of decreasing size (500 µm, 200 µm and 80 µm) and each sieve content was then collected onto pre-weighed and pre-combusted GF/F filters. A second net sample was kept in a large volume of water until the return to the laboratory. Different groups and taxa were then sorted directly under a dissecting microscope: Acartia tonsa (copepod Calanoida), harpacticoid copepods, cirriped nauplii, mysids, Monocorophium insidiosum (amphipod), isopods, and gastropod larvae. About 300 to 500 individuals from each taxa were sorted, placed together on pre-combusted GF/F filters, rinsed with ultrapure water and kept in hermetically sealed boxes. Cnidarians and ctenophores were collected with a dip net, rinsed with ultrapure water, measured (nearest mm) and individually placed in polyethylene bags. To minimize bias in C and N isotope analysis related to the presence of prey in their digestive tract, zooplankton and gelatinous organisms were conditioned only three hours after sampling, enabling gut clearing. Juveniles of teleostean fish (Mugil sp.) were caught with a landing net and measured (nearest mm). Dorsal white muscle samples without skin were individually placed in polyethylene bags. All samples were stored in a freezer at -20 • C before being freeze-dried for 24 h, then stored in a dry place. Numbers of samples analyzed are given in Table 2.

Stable isotope analysis

Before stable isotope analysis, all freeze-dried samples were ground using a mortar and pestle and sieved on a steel sieve (60 µm) to obtain a homogeneous fine powder. Several compartments (POM, SOM, macrophytes, crustaceans, micro-and mesozooplankton, gastropod larvae), susceptible to contain carbonates, were split into two subsamples. One sub-sample was acidified to remove carbonates, as they are enriched in 13 C [START_REF] Lorrain | Decarbonation and preservation method for the analysis of organic C and N contents and stable isotope ratios of low-carbonated suspended particulate material[END_REF] and was used for measurements of δ 13 C values. The other subsample was used for measurements of δ 15 N values of the organisms without any prior treatment. SOM decarbonation was done with 1% HCl solution abundantly poured onto the sample until cessation of bubbling. POM samples were acidified by three successive rinses of filters with 1% HCl and then rinsed thoroughly with ultrapure water. For macrophytes, crustaceans, zooplankton and gastropod larvae, acidification was carried out by exposure to HCl fumes for four hours. All samples were then rinsed three times with ultrapure water and dried in an oven (60 • C). Non-acidified samples were placed in tin microcapsules, while decarbonated samples were conditioned in silver microcapsules. The mass of material analyzed was about 0.5 mg for crustaceans, gelatinous organisms and teleosteans, and 1 mg for POM and SOM. Stable isotope analyses (C and N) were carried out at the stable isotope facility of the LIttoral ENvironnement et Sociétés (LIENSs) stable isotope facility (CNRS -University of La Rochelle), using an elemental analyzer (Flash EA 1112, Thermo Scientific) coupled with an isotope ratio mass spectrometer (Delta V Advantage with a Conflo IV interface, Thermo Scientific).

Isotope composition is expressed in the δ notation ( ) as deviations from international standards (Vienna Pee Dee belemnite for carbon and atmospheric N 2 for nitrogen, [START_REF] Peterson | Stable isotopes in ecosystem studies[END_REF]: δX = [(Rsample/Rstandard) -1] ×10 3 , with X is 13 C or 15 N and R the isotopic ratio 13 C/ 12 C or 15 N/ 14 N. Carbon/Nitrogen (C/N) atomic ratios were also computed. Calibration was done using reference materials (USGS-24, IAEA-CH6, IAEA-600, USGS-61, USGS-62 for carbon; IAEA-N2, IAEA-NO-3 and IAEA-600, USGS-61, USGS-62 for nitrogen). Analytical precision based on the analyses of USGS-61 and USGS-62 used as laboratory internal standard was <0.10 and <0.15 for carbon and nitrogen, respectively.

As lipids are depleted in 13 C compared to the other biochemical components, a high lipid content may have a significant effect on δ 13 C values of organisms [START_REF] De Niro | Mechanism of carbon isotope fractionation associated with lipid synthesis[END_REF][START_REF] Logan | Lipid corrections in carbon and nitrogen stable isotope analyses: comparison of chemical extraction and modelling methods[END_REF] and a correction should be applied for zooplanktonic organisms presenting a C/N ratio >5 [START_REF] Smyntek | A standard protocol for stable isotope analysis of zooplankton in aquatic food web research using mass balance correction models[END_REF]. Mnemiopsis leidyi δ 13 C values were thus corrected for lipids according = 1.3 ).

Data analysis

Differences in mean stable isotopic values and C/N ratio among food web components or organism size classes were tested using nonparametric Kruskall-Wallis ANOVA, followed by post-hoc Dunn tests of comparison of means. A hierarchical clustering analysis, using Euclidean distance and Ward criterion, enabled the grouping of species of similar isotope composition. The relationships between the logarithm of individual size and δ 15 N value in M. leidyi and Aurelia sp. were tested by Pearson's linear correlation. The respective contributions of the different sources of organic matter ingested by M. leidyi and Aurelia sp. were determined with Bayesian mixing models using the SIAR (Stable Isotope Analysis in R) package [START_REF] Parnell | Source partitioning using stable isotopes: Coping with too much variation[END_REF] although these models are called into question below [START_REF] Stock | Analyzing mixing systems using a new generation of Bayesian tracer mixing models[END_REF]. SIAR models estimate probability distributions of the contribution of different sources to consumers. The Bayesian approach of mixing models takes into account the variability in stable isotope values of sources and provides a probability distribution associated with contribution values, but it is sensitive to the choice of the trophic enrichment factor (TEF) applied [START_REF] Phillips | Best practices for use of stable isotope mixing models in food-web studies[END_REF][START_REF] Bond | Recent Bayesian stable-isotope mixing models are highly sensitive to variation in discrimination factors[END_REF]. The TEF values used for the gelatinous organisms were selected following [START_REF] Vansteenbrugge | The Non-Indigenous Ctenophore Mnemiopsis Leidyi in the Southern North Sea: Ecological and Socio-Economic Effects Related To Its Trophic Position and the Current Distribution of Gelatinous Zooplankton[END_REF]: ∆δ 13 C = 1.6 ± 0.6 and ∆δ 15 N = 3.5 ± 0.6 , for M. leidyi. SIAR models were run with two prey categories (Plankton and Benthos) according to the results of a previous study on the gastrovascular cavity contents of M. leidyi (Marchessaux et al., submitted). However, as different authors testified from predation on ctenophores by jellyfishes [START_REF] Purcell | Predation by the scyphomedusan Chrysaora quinquecirrha on Mnemiopsis leidyi ctenophores[END_REF][START_REF] Breitburg | Predator-mediated landscape structure: seasonal patterns of spatial expansion and prey control by Chrysaora quinquecirrha and Mnemiopsis leidyi[END_REF], three prey categories (Plankton, Benthos and Mnemiopsis) were also considered as potential prey for Aurelia sp. in SIAR models.

Results

Stable isotope composition of food web components

Stable isotope values of planktonic (POM, micro-and mesozooplankton, ctenophores and cnidarians) and benthic (SOM, macroalgae, seagrasses, crustaceans and juvenile fish) components analyzed in Berre Lagoon ranged from -27.1 to -9.5 for δ 13 C values and from 9.2 to 15.1 for δ 15 N values (Table 1). At the base of the food web, the wide range of δ 13 C values (∆δ 13 C = 14.6 ) clearly distinguished the seagrass Zostera noltii (ZOS) from the other carbon sources (Chaetomorpha linum (CHA), POM and SOM), which presented much lower δ 13 C values (from -25.9 to -22.3 ). All these organic matter (OM) sources presented low δ 15 N values (from 9.2 to 10.1 ). The range of δ 15 N (∆δ 15 N = 5.9 ) between the base of the food web and Aurelia sp. (AUR: 15.1 ) suggested the existence of two trophic levels between them.

The hierarchical clustering separated two groups of potential prey for gelatinous organisms, which presented significantly different carbon and nitrogen compositions (Fig. 2): benthic prey (δ 13 C = -19.8 ± 0.9 and δ 15 N = 12.1 ± 0.3 ) and planktonic prey (δ 13 C = -26.1 ± 1.0 and δ 15 N = 11.5 ± 0.4 ). The group of benthic organisms (hereafter Benthos) included harpacticoid copepods (HAR), isopods (ISO), and the amphipod Monocorophium insidiosum (COR), and displayed significantly higher δ 13 C and δ 15 N values (ANOVA, F = 317.30, p < 0.001 and F = 18.53, p < 0.001 for C and N respectively) than the group of planktonic prey (hereafter Plankton) composed by microzooplankton 80-200 µm (Z80), mesozooplankton >200 µm (Z200), Acartia tonsa (ACA), cirriped nauplii (CIR) and gastropod larvae (GAS). Among planktonic organisms, Z80 and GAS presented close C and N isotope composition with the lowest δ 15 N values (<11.5 ), while ACA, CIR and Z200 exhibited higher, while not significantly, δ 15 N values (>11.5 ). Juveniles of the teleost fish Mugil sp. (MUG) had a similar isotope composition to the group of benthic prey, but due to their size (>5 cm) were not potential prey for gelatinous organisms and were thus excluded from the Benthos prey group. Mysid shrimps (MYS), abundant in Zostera seagrass meadows, displayed a δ 13 C value (-9.5 ± 0.2 ) intermediate between those of Plankton and Benthos, and a significantly higher δ 15 N value (Table 1, Fig. 2), but were never observed in M. leidyi gastrovascular cavity (GM pers. obs.). M. leidyi (MNE) was positioned close to the Plankton group with slightly higher δ 13 C and δ 15 N values, while Aurelia sp. (AUR) presented the highest δ 15 N values (δ 15 N = 15.1 ± 0.9 ) observed among the organisms analyzed (Fig. 2).

Size-and species-related variations

Three size classes were distinguished among MNE individuals collected, which were all of small size: 10 mm, 15 mm and 20 mm, which significantly differed in their stable isotope composition (Fig. 3): individuals from the 20 mm size class were significantly more enriched in 13 C than the two smaller classes (-23.2 vs. <-25.0 , KW-ANOVA, H = 9.44, p = 0.009) and individuals from the smallest 10 mm size class were significantly more depleted in 15 N than the two largest ones (10.2 vs. >12.0 , H = 8.84, p = 0.012). A significant positive linear correlation was observed between δ 15 N value and the logarithm of size in the two species (r = 0.70, p = 0.001, N = 18 and r = 0.77, p = 0.015, N = 9, for M. leidyi and Aurelia sp., respectively).

Two categories (Plankton and Benthos) with statistically different stable isotope composition were identified as potential prey for both M. leidyi and Aurelia sp. and used in SIAR models. Results testified to a zooplanktonic based diet for the two gelatinous organisms, but also to the assimilation of a small proportion of benthic organisms (Fig. 4). A mean proportion of 97.5% (93.6-100%, 95% Bayesian credible interval) of Plankton was found in MNE diet, with a low but significant proportion of Benthos (mean 2.5%; 0.3-6.4%) (Fig. 4a). SIAR model with two prey categories indicated that AUR mainly preyed on plankton (79.2%; 68.6-89.9%) but ingested also a small quantity of benthic organisms or some prey influenced by benthic OM (20.8%; 10.1-31.4%) (Fig. 4b). When MNE was included in the SIAR model as a third prey category for AUR, the mean proportions of assimilated prey in AUR tissues were 47.9% (26.1-71.8%, 95% Bayesian credible interval) for Plankton, 12.1% (10.2-66.5%) for Benthos and 39.9% (18.8-22.0%) for MNE.

Discussion

Isotopic composition of jellyfish food webs in Berre Lagoon

The carbon and nitrogen isotope composition of M. leidyi and Aurelia sp. indicated that their food resources were mainly based on planktonic POM, and in a lower proportion on benthic SOM and macroalgae, while the seagrass Zostera noltii was not used as a food source directly (Fig. 2). The high δ 13 C value of Z. noltii in Berre Lagoon (-9.5 ) is similar to that (-9.7 ) reported for this seagrass by [START_REF] Carlier | Spatial heterogeneity in the food web of a heavily modified Mediterranean coastal lagoon: stable isotope evidence[END_REF] in Salses-Leucate, a lagoon located in the southwestern part of the French Mediterranean coast. POM includes phytoplankton, but also microzooplankton, bacteria, detritus, and fecal pellets, among other elements [START_REF] Lienart | Dynamics of particulate organic matter composition in coastal systems: A spatio-temporal study at multi-systems scale[END_REF]. In coastal waters POM presents generally lower δ 13 C and higher δ 15 N values than phytoplankton, mainly due to the influence of terrestrial detritus in its composition, as observed in the Gulf of Lions [START_REF] Harmelin-Vivien | Comparison of C and N stable isotope ratios between surface particulate organic matter and microphytoplankton in the Gulf of Lions (NW Mediterranean)[END_REF]. Berre Lagoon is exposed to high freshwater inputs and thus presents a low δ 13 C value (-25.9 ), close to those of the Rhône River POM (-27.4 to -26.6 ;[START_REF] Harmelin-Vivien | Seasonal variation in stable C and N isotope ratios of the Rhone River inputs to the Mediterranean Sea (2004-2005)[END_REF]. Low δ 13 C (< -24 ) and relatively high δ 15 N (> 9 ) values for POM, as observed in Berre Lagoon, were reported in other highly anthropized French Mediterranean lagoons subject to river inputs or with restricted connections with restricted exchanges with the sea, such as Mauguio [START_REF] Vizzini | Spatial variability of stable carbon and nitrogen isotope ratios in a Mediterranean coastal lagoon[END_REF][START_REF] Dierking | Fingerprints of lagoonal life: migration of the marine flatfish Solea solea assessed by stable isotopes and otolith microchemistry[END_REF] and Canet [START_REF] Carlier | Spatial heterogeneity in the food web of a heavily modified Mediterranean coastal lagoon: stable isotope evidence[END_REF]. In other French and Italian lagoons subject to lower river inputs or benefiting from more extensive connections with the sea, slightly higher δ 13 C (-23 to -18 ) and lower δ 15 N (2 to 6 ) values were observed for POM [START_REF] Carlier | A seasonal survey of the food web in the Lapalme Lagoon (northwestern Mediterranean) assessed by carbon and nitrogen stable isotope analysis[END_REF][START_REF] Vizzini | The fate of organic matter sources in coastal environments: a comparison of three Mediterranean lagoons[END_REF]. In all the studied Mediterranean lagoons, SOM had higher (or similar) δ 13 C values than POM, as observed in Berre Lagoon, highlighting a similar origin, while varying δ 15 N values were more related to the degree of human impact. SOM was composed of settled POM, microphyto-and meio-benthos, and a wide range of detritus (detrital macroalgae, seagrasses and terrestrial plants, fecal pellets, etc.) [START_REF] Cresson | Spatiotemporal variation of suspended and sedimentary organic matter quality in the Bay of Marseilles (NW Mediterranean) assessed by biochemical and isotopic analyses[END_REF]. In Berre lagoon, the higher δ 13 C value of SOM was also probably due to the influence of Z. noltii detritus.

Two groups of potential prey items for ctenophores and jellyfishes with statistically different stable isotope values were individualized in Berre Lagoon. The pelagic group (Z80, Z200, the copepod A. tonsa, cirriped nauplii and gastropod larvae) included planktonic organisms relying on POM and phytoplankton as carbon sources, as observed by [START_REF] Vansteenbrugge | Trophic ecology of Mnemiopsis leidyi in the southern North Sea: a biomarker approach[END_REF] in the North Sea. Microzooplankton (Z80) and gastropod larvae displayed the lowest δ 15 N values as they mainly fed on picoand nanophytoplankton, heterotrophic protozoans and bacteria [START_REF] Calbet | Phytoplankton growth, microzooplankton grazing, and carbon cycling in marine systems[END_REF][START_REF] Vargas | Feeding by larvae of intertidal invertebrates: assessing their position in pelagic food webs[END_REF]. Mesozooplankton (Z200), cirriped larvae and the invasive copepod A. tonsa exhibited slightly higher δ 15 N values, while not statistically different (Table 2), as they fed on larger cells (diatoms, chlorophytes, dinoflagellates, tintinnids and rotifers) [START_REF] Gaudy | Copepod nutrition in a Mediterranean lagoon in terms of particle concentration and temperature[END_REF][START_REF] Saiz | Feeding rates and prey: predator size ratios of the nauplii and adult females of the marine cyclopoid copepod Oithona davisae[END_REF].

The second group of potential prey items included benthic organisms such as harpacticoid copepods, the amphipod Monocorophium insidiosum and isopods, which presented significantly higher δ 13 C values than planktonic prey, and slightly higher δ 15 N values. [START_REF] Vizzini | Seasonal variations in the stable carbon and nitrogen isotope ratios ( 13 C/ 12 C and 15 N/ 14 N) of primary producers and consumers in a western Mediterranean coastal lagoon[END_REF] also recorded higher δ 13 C values in harpacticoid copepods and Monocorophium amphipods than in planktonic calanoid copepods in Sabaudia Lagoon. These high δ 13 C values were in accordance with a consumption of benthic OM (SOM, benthic diatoms and macrophytes detritus) by these compartments [START_REF] Nair | Life cycle of Corophium insidiosum (Crustacea, Amphipoda) in laboratory culture[END_REF][START_REF] Arrontes | Diet, food preference and digestive efficiency in intertidal isopods inhabiting macroalgae[END_REF]. As observed in Berre Lagoon, rather high δ 13 C and δ 15 N values were recorded for mysids in different Mediterranean lagoons and in the North Sea coastal areas [START_REF] Persic | Stable carbon and nitrogen isotope quantitative structural assessment of dominant species from the Vaccarès Lagoon trophic web (Camargue Biosphere Reserve, France)[END_REF][START_REF] Vizzini | The fate of organic matter sources in coastal environments: a comparison of three Mediterranean lagoons[END_REF][START_REF] De Clippele | Assessing the Trophic Position of the Invasive Ctenophore Mnemiopsis Leidyi in the Southern North Sea[END_REF]. The positioning of these crustaceans at a higher trophic level than planktonic and benthic smaller crustaceans was in accordance with their predatory feeding behavior [START_REF] Rudstam | Diel vertical migration and feeding patterns of Mysis mixta (Crustacea, Mysidacea) in the Baltic Sea[END_REF].

The δ 13 C values (-26 to -22 ) observed for M. leidyi in Berre lagoon were similar to those measured by [START_REF] Vansteenbrugge | Trophic ecology of Mnemiopsis leidyi in the southern North Sea: a biomarker approach[END_REF] in the port of Dunkirk (-24 to -22 ), but lower than those generally recorded in the coastal and estuarine sites of the North Sea (-22 to -16 ) [START_REF] Hamer | The invasive ctenophore Mnemiopsis leidyi: a threat to fish recruitment in the North Sea?[END_REF][START_REF] De Clippele | Assessing the Trophic Position of the Invasive Ctenophore Mnemiopsis Leidyi in the Southern North Sea[END_REF][START_REF] Vansteenbrugge | Trophic ecology of Mnemiopsis leidyi in the southern North Sea: a biomarker approach[END_REF][START_REF] Van Walraven | Temporally varying isotopic niche overlap of the invasive ctenophore Mnemiopsis leidyi with other zooplanktivores in the Dutch Wadden Sea[END_REF].

The δ 15 N values for this ctenophore (9 to 13 ) in Berre Lagoon were lower than those reported by these authors in the North Sea (12 to 21 ), probably owing to the small size of the individuals collected. However, similar C and N isotopic composition was recorded for Aurelia sp. in Berre Lagoon (δ 13 C = -23 to -23.5 ; δ 15 N = 13 to 16 ) and the North Sea (δ 13 C = -24 to -18.5 ; δ 15 N = 11 to 15 ) for individuals of similar size [START_REF] Javidpour | Temporal dietary shift in jellyfish revealed by stable isotope analysis[END_REF][START_REF] Van Walraven | Temporally varying isotopic niche overlap of the invasive ctenophore Mnemiopsis leidyi with other zooplanktivores in the Dutch Wadden Sea[END_REF].

Variation in ctenophore and jellyfish diet: planktonic vs. benthic prey

The two gelatinous organisms investigated in Berre Lagoon seemed to be positioned at different trophic levels, but this was mainly due to the considerable difference of size between the (B. Hunt, pers. comm.), close to those observed here for Aurelia sp. of similar size (13.6 to 16.2 ; 80-160 mm), positioning the two species at a similar trophic level for similar sizes. Thus, the difference and narrow range of C and N isotope composition recorded for these two jellyfishes in this study most likely came from the low number and narrow size range of the individuals analyzed. In the future, it would be important to consider a wide range of ctenophore and cnidarian size classes in order to determine with precision their role in the food web functioning of Berre Lagoon.

Significant differences in both isotope compositions were observed among the three size classes of M. leidyi, with the smallest size class presenting the lowest δ 13 C and δ 15 N values (Fig. 3). This suggested that 10 mm-sized M. leidyi probably consumed phytoplankton (POM) and small microzooplanktonic organisms. Small ctenophores (10 mm) are reported to be omnivorous, consuming small prey items such as tintinnids, protists and copepod nauplii [START_REF] Sullivan | Diet of the larval ctenophore Mnemiopsis leidyi A. Agassiz (Ctenophora, Lobata)[END_REF][START_REF] Rapoza | Life-stage dependent, in situ dietary patterns of the lobate ctenophore Mnemiopsis leidyi Agassiz 1865[END_REF], and to be positioned at a low trophic level in the food webs. The largest M. leidyi size class analyzed in the present study was only 20 mm long, but it presented significantly higher δ 13 C and δ 15 N values than the smaller ones, indicating a shift in their diet toward larger prey (micro-and mesozooplankton) and benthic crustaceans, as suggested by [START_REF] Vansteenbrugge | Trophic ecology of Mnemiopsis leidyi in the southern North Sea: a biomarker approach[END_REF]. This dietary pattern was confirmed by the results of the analysis of the gastrovascular cavity content of M. leidyi in Berre Lagoon (Marchessaux et al., submitted), indicating that this ctenophore ingested mainly mollusk larvae (34%), copepods (29%), barnacle larvae (13%) copepod nauplii (11%) and the benthic amphipod M. insidiosum.

Although exploratory in this study and the subject of criticism, the mixing models SIAR allowed estimation of the contributions of the different food sources assimilated by the two gelatinous organisms. SIAR results indicated that planktonic prey dominated the diet of both species, particularly the copepod A. tonsa, cirriped nauplii, and gastropod larvae, clustered in the Plankton group. The isotope composition of the ctenophore close to that of its prey was mainly explained by a very rapid assimilation rate, particularly for small individuals such as those analyzed in this study [START_REF] Kremer | Growth dynamics of a ctenophore (Mnemiopsis) in relation to variable food supply. II. Carbon budgets and growth model[END_REF][START_REF] Reeve | Growth dynamics of a ctenophore (Mnemiopsis) in relation to variable food supply, I. Carbon biomass, feeding, egg production, growth and assimilation efficiency[END_REF][START_REF] Bumann | The ctenophore Mnemiopsis leidyi has a flow-through system for digestion with three consecutive phases of extracellular digestion[END_REF]. However, SIAR results indicated that both M. leidyi and Aurelia sp. also assimilated benthic prey in addition to pelagic prey, even if in small proportions (from 2.5% to 21%; Fig. 4. The predation on benthic organisms by these two gelatinous organisms was confirmed by the analysis of their gastrovascular cavity content, within which harpacticoid copepods and small Monocorophium insidiosum were often recorded [START_REF] Marchessaux | Etude de L'Anthroposystème Emblématique de L'étang de Berre : Approches écosystémique Et Sociologique de L'Impact Du Cténaire Invasif Mnemiopsis Leidyi[END_REF]Marchessaux et al., submitted)). The predation on benthic organisms by pelagic ctenophores and cnidarians was due to specific behavior patterns of both predator and prey. M. leidyi was observed swimming close to the bottom or among seagrasses when the wind was strong and the surface water rough, or at night (GM pers. obs.). Such behavior favored its predation on the small crustaceans (i.e. amphipods, harparcticoids) from the sediment, as observed by [START_REF] Pitt | Predation by jellyfish on large and emergent zooplankton: Implications for benthic-pelagic coupling[END_REF], or living in the Zostera canopy. In addition, juvenile M. insidiosum, which were among the most common amphipods found in the Mnemiopsis gastrovascular cavity, are known for being very active outside their tubes in contrast to the adults [START_REF] Fricke | Insights into ecological and reproductive aspects of two cryptogenic peracarid crustaceans of the Argentinian coast[END_REF]. The predation on benthic crustaceans, such as the Malacostraca Lucifer and night emergent copepods, by the scyphozoan jellyfish Catostylus mosaicus in Australia was demonstrated by [START_REF] Pitt | Predation by jellyfish on large and emergent zooplankton: Implications for benthic-pelagic coupling[END_REF], providing evidence of benthic-pelagic coupling by jellyfishes in coastal areas. Evidence of assimilation of benthic OM by jellyfishes was also provided through the use of sulfur stable isotope (δ 34 S) in Aurelia aurita in the North Sea [START_REF] Javidpour | Temporal dietary shift in jellyfish revealed by stable isotope analysis[END_REF] and of fatty acids in A. aurita and two other jellyfish species in the Yellow Sea [START_REF] Cui | Potential dietary influence on the stable isotopes and fatty acid compositions of jellyfishes in the Yellow Sea[END_REF]. Thus, pelagic cnidarians and ctenophores could play a major and underestimated role in benthic-pelagic coupling of food webs in coastal areas.

Ctenophore-jellyfish relationships: competition and/or predation?

Particulate Organic Matter (POM) was the main OM source for zooplankton and thus the major base line of the planktonic food web of Berre Lagoon. C and N stable isotope composition and gastrovascular cavity content indicated that M. leidyi and Aurelia sp. mainly relied on the same planktonic prey for food (mollusk veliger larvae, particularly gastropod larvae, copepods, mainly A. tonsa, cirriped nauplii, the cladoreran Pleopis polyphaemoides, and polychete larvae) (Marchessaux et al., submitted), being thus in potential competition when these resources decreased. [START_REF] Hamer | The invasive ctenophore Mnemiopsis leidyi: a threat to fish recruitment in the North Sea?[END_REF] also reported potential competition between the invasive Mnemiopsis and the native ctenophore Pleurobrachia pileus.

But predation relationships among pelagic ctenophores and cnidarians have been reported by different authors [START_REF] Kreps | Escape of the ctenophore Mnemiopsis leidyi from the scyphomedusa predator Chrysaora quinquecirrha[END_REF][START_REF] Purcell | Predation by the scyphomedusan Chrysaora quinquecirrha on Mnemiopsis leidyi ctenophores[END_REF][START_REF] Hosia | Interactions between native and alien ctenophores: Beroe gracilis and Mnemiopsis leidyi in Gullmarsfjorden[END_REF][START_REF] Tilves | Predation by the scyphozoan Pelagia noctiluca on Mnemiopsis leidyi ctenophores in the NW Mediterranean Sea[END_REF][START_REF] Tilves | Natural diet and predation impacts of Pelagia noctiluca on fish eggs and larvae in the NW Mediterranean[END_REF][START_REF] Breitburg | Predator-mediated landscape structure: seasonal patterns of spatial expansion and prey control by Chrysaora quinquecirrha and Mnemiopsis leidyi[END_REF][START_REF] Breitburg | Predator-mediated landscape structure: seasonal patterns of spatial expansion and prey control by Chrysaora quinquecirrha and Mnemiopsis leidyi[END_REF][START_REF] Frost | Distribution and trophic links of gelatinous zooplankton on Dogger Bank, North Sea[END_REF] confirmed Beroe sp. (an indigenous ctenophore) as a predator of M. leidyi, based on isotope composition. In Berre Lagoon, the cnidarian Aurelia sp. was positioned at a higher trophic level than the ctenophore, probably owing to the small size of the M. leidyi individuals analyzed, and SIAR results suggested that small Mnemiopsis could also represent potential prey for large Aurelia sp. Predation on Aurelia planula and ephyrae by M. leidyi was observed in the Baltic Sea [START_REF] Javidpour | Annual assessment of the predation of Mnemiopsis leidyi in a new invaded environment, the Kiel Fjord (Western Baltic Sea): a matter of concern?[END_REF][START_REF] Riisgård | High abundance of the jellyfish Aurelia aurita excludes the invasive ctenophore Mnemiopsis leidyi to establish in a shallow cove (Kertinge Nor, Denmark)[END_REF]. It might be possible that the opposite is also true and that Aurelia adults were preying on the ctenophore eggs, larvae, and juveniles, as suggested here. However, this hypothesis needs to be thoroughly tested in situ and with experiments.

Conclusion

While exploratory, the present study highlighted a potential important predation role of the invasive ctenophore Mnemiopsis leidyi and the native moon jellyfish Aurelia sp. on the same zooplanktonic organisms in Berre Lagoon. Our results evidenced that both M. leidyi and Aurelia sp. preyed on and also assimilated benthic organisms (harpacticoid copepods, amphipods), while in lower proportions, which confirmed the results of the gastrovascular cavity content analysis (Marchessaux et al., submitted). The shallow mean depth of Berre Lagoon probably favored two pathways of OM transfer (i.e. from both pelagic and benthic communities) in which pelagic ctenophores and cnidarians played an underestimated role. In the future, a more complete study of the functioning of the whole food web, up to fishes, should be carried out in this lagoon, including in particular an assessment of the seasonal and spatial variations of the C and N isotope composition of both the organisms and the organic matter sources.

First

  overview on trophic relationships of the invasive ctenophore Mnemiopsis leidyi in a Mediterranean coastal lagoon (Berre Lagoon, France): benthic-pelagic coupling evidenced by carbon and nitrogen stable isotope composition Guillaume Marchessaux, Mireille Harmelin-Vivien, M. Ourgaud, D. Bănaru, Loic Guilloux, Bruno Belloni, Benoit Lebreton, Gaël Guillou, Delphine Thibault To cite this version: Guillaume Marchessaux, Mireille Harmelin-Vivien, M. Ourgaud, D. Bănaru, Loic Guilloux, et al.. First overview on trophic relationships of the invasive ctenophore Mnemiopsis leidyi in a Mediterranean coastal lagoon (Berre Lagoon, France): benthic-pelagic coupling evidenced by carbon and nitrogen stable isotope composition. Regional Studies in Marine Science, Elsevier, 2021, 41, pp.101570. 10.1016/j.rsma.2020.101570. hal-03362598

Fig. 1 .

 1 Fig. 1. Study area and sampling station. SP: sampling station. Black triangles: freshwater inputs (Touloubre, Durance (EDF powerplant), Arc) and white triangle: connection to the sea (Caronte Canal, south-west). Med. Sea: Mediterranean Sea.

Fig. 2 .

 2 Fig. 2. Biplot of mean (± SD) δ 15 N and δ 13 C ( ) isotope values in the jellyfish food web of Berre Lagoon in April 2016. ACA: Acartia tonsa; AUR: Aurelia sp. (80-160 mm) ; CHA: Chaetomorpha linum; CIR: cirriped nauplii; COR: juveniles of Monocorophium insidiosum; GAS: gastropod larvae; HAR: harpacticoids; ISO: isopods; MNE: Mnemiopsis leidyi (10-20 mm); MUG: juveniles of Mugil sp.; MYS: mysids; POM: particulate organic matter; SOM: surface sediment; Z80: zooplankton 80-200 µm; Z200: zooplankton 200-500 µm; ZOS: Zostera noltii.

  to the equation given by[START_REF] Post | Getting to the fat of the matter: models, methods and assumptions for dealing with lipids in stable isotope analyses[END_REF], as done by D'Ambra et al.(2014, 2018) for Aurelia sp.. Mean change in δ 13 C values (M. leidyi ∆δ 13 C = 1.9 ) due to lipid correction was close to the value given by[START_REF] Logan | Lipid corrections in carbon and nitrogen stable isotope analyses: comparison of chemical extraction and modelling methods[END_REF] for marine invertebrates (∆δ 13 C

Fig. 3 .

 3 Fig. 3. Boxplot of mean (± SD) isotope values of Mnemiopsis leidyi, (A) δ 13 C ( ) and (B) δ 15 N ( ), as a function of ctenophore size (mm). Number of individuals: 10 mm (N = 4), 15 mm (N = 10), 20 mm (N = 4). Horizontal bars represent median values and red crosses mean values.

Fig. 4 .

 4 Fig. 4. Relative contributions of planktonic (Plankton) and benthic (Benthos) prey to the diets of Mnemiopsis leidyi and Aurelia sp. from SIAR models. Points represent the mode of the distribution and boxes the 50%, 75% and 95% Bayesian credibility intervals from light to dark gray.. individuals analyzed in the two species. The small-sized (10 to 20 mm) M. leidyi analyzed presented a lower δ 15 N value (12.1 ) than the large-sized (80 to 160 mm) Aurelia sp. (15.1 ). However, a significant positive correlation between δ 15 N value and size was observed in the two species in this lagoon, as recorded by Javidpour et al. (2016) for Aurelia aurita in the Baltic Sea, and Van Walraven et al. (2018) for M. leidyi in the North Sea. But Vansteenbrugge et al. (2016) observed unexpectedly higher δ 15 N values in small (< 10 mm) than in large (11 to 55 mm) M. leidyi individuals in the Westerschelde estuary. Higher δ 15 N values (16.1 to 16.5 ) have been observed for large M. leidyi (70-120 mm) in Berre Lagoon in a previous study(B. Hunt, pers. comm.), close to those observed here for Aurelia sp. of similar size (13.6 to 16.2 ; 80-160 mm), positioning the two species at a similar trophic level for similar sizes. Thus, the difference and narrow range of C and N isotope composition recorded for these two jellyfishes in this study most likely came from the low number and narrow size range of the individuals analyzed. In the future, it would be important to consider a wide range of ctenophore and cnidarian size classes in order to determine with precision their role in the food web functioning of Berre Lagoon.Significant differences in both isotope compositions were observed among the three size classes of M. leidyi, with the smallest

Table 1

 1 Variability of environmental parameters in the Berre Lagoon with mean values during the present study, and minimum-maximum values between 2010 and 2017[START_REF] Marchessaux | Refugia area for the ctenophore Mnemiopsis leidyi A. Agassiz 1865 in the Berre Lagoon (southeast France): The key to its persistence[END_REF].

	Water parameters	Mean ± SD (April 2016) Min-Max (2010-2017)
	Temperature ( • C)	16.5 ± 0.3	0.3 -28.2
	Salinity	17.4 ± 0.8	14.6 -33.7
	Chlorophyll a (µg L -1 )	1.1 ± 0.8	0.2 -37.1
	Oxygen (µmol L -1 )	245.1 ± 2.4	1.7 -385.6

Table 2

 2 Mean (± standard deviation) δ 15 N and δ 13 C (‰) values of the different compartments analyzed in the Berre Lagoon food webs (NW Mediterranean). The letters indicate whether the differences are significant between the compartments (Kruskal-Wallis ANOVA. p < 0.05). N = number of samples.

	Compartments	Codes	N	δ 13 C				δ 15 N				C / N
	Organic matter sources
	POM	POM	4	-25.9	±	0.1	a,b	9.9	±	0.4	a,b	3.9	±	1.0	a,b
	SOM	SOM	4	-22.3	±	0.0	d,e	10.1	±	0.1	a,b	8.3	±	0.1	a,b
	Zostera noltii	ZOS	3	-9.5	±	0.2	e	9.4	±	0.2	a	24.4	±	0.7	b
	Chaetomorpha linum	CHA	3	-24.1	±	0.1	b,c,d,e	9.2	±	0.1	a	8.1	±	0.2	a,b
	Zooplankton
	Zooplankton 80-200 µm	Z80	15	-25.5	±	1.2	a,b	11.1	±	0.1	a,b	2.8	±	0.8	a
	Zooplankton > 200 µm	Z200	15	-26.3	±	0.4	a	11.8	±	0.2	c	3.7	±	0.6	a,b
	Acartia tonsa	ACA	7	-27.1	±	0.8	a	11.7	±	0.2	b,c	4.0	±	0.2	a,b
	Cirripeds nauplii	CIR	4	-25.5	±	0.4	a,b,c	11.7	±	0.1	b,c	4.7	±	0.7	a,b
	Gastropod larvae	GAS	1	-25.6				11.3				4.0
	Gelatinous zooplankton
	Aurelia sp. (80-160 mm)	AUR	9	-23.2	±	0.1	c,d,e	15.1	±	0.9	e	4.2	±	0.1	a,b
	Mnemiopsis leidyi (10-20 mm)	MNE	18	-24.8	±	1.0	b,c,d	12.0	±	1.2	c,d	5.4	±	2.2	a,b
	Benthic invertebrates
	Monocorophium insidiosum	COR	4	-20.7	±	0.4	e	11.9	±	0.2	c,d	4.7	±	0.3	a,b
	Harparcticoids	HAR	3	-19.6	±	0.4	e	12.4	±	0.2	c,d,e	4.3	±	0.0	a,b
	Isopods	ISO	3	-18.9	±	0.4	e	12.0	±	0.1	c,d	4.1	±	1.0	a,b
	Mysids	MYS	3	-22.9	±	0.3	c,d,e	13.8	±	0.2	d,e	3.6	±	0.1	a,b
	Teleosteans
	Mugil sp. Juveniles	MUG	6	-21.5	±	0.9	e	12.2	±	0.9	c,d	3.8	±	0.2	a,b
	Kruskal-Wallis ANOVA			H = 82.8 p < 0.001		H = 73 .1 p < 0.001		H = 31.3 p < 0.005

Acknowledgment

We would like to thank Robert Arfi (MIO), Raphaël Grisel (GIPREB), Adrien Goujard (GIS Posidonie), Dorian Guillemain (MIO) and Elodie Rouanet (GIS Posidonie) for their help during in situ sampling. Thanks to Maria Jesus Ferrara-Guerrero (Universidad Nacional Autónoma de México), Maria-Elena Huidobro-Salas (Universidad Nacional Autónoma de México), José Luis Gama-Flores (Universidad Nacional Autónoma de México) and Marc Pagano (MIO) for their help in sample conditioning. We also thank the Microscopie et Imagerie Marine (MIM) and Service à la Mer (SAM) facilities in MIO. The project leading to this publication has received funding from the European FEDER Fund under project 1166-39417. Guillaume Marchessaux was supported by a PhD fellowship from the French Ministry of Higher Education and Research from 2015 to 2018. Thanks to Michael Paul, a native English speaker, for proof-reading the text.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.


Chunxia Qin

Puxun Tu

Xiaojun Chen
email: xiaojunchen@sjtu.edu.cn.

Jocelyne Troccaz

A novel registration-based algorithm for prostate segmentation via the combination of SSM and CNN

Keywords: registration-based segmentation, statistical shape mode, probability map, boundary distance map i

Purpose: Precise determination of target is an essential procedure in prostate interventions, such as prostate biopsy, lesion detection, and targeted therapy. However, the prostate delineation may be tough in some cases due to tissue ambiguity or lack of partial anatomical boundary. In this study, we proposed a novel supervised registrationbased algorithm for precise prostate segmentation, which combine the convolutional neural network (CNN) with a statistical shape model (SSM). Methods: The proposed network mainly consists of two branches. One called SSM-Net branch was exploited to predict the shape transform matrix, shape control parameters, and shape fine-tuning vector, for the generation of the prostate boundary. Furtherly, according to the inferred boundary, a normalized distance map was calculated as the output of SSM-Net. Another branch named ResU-Net was employed to predict a probability label map from the input images at the same time. Integrating the output of these two branches, the optimal weighted sum of the distance map and the probability map was regarded as the prostate segmentation. Results: Two public datasets PROMISE12 and NCI-ISBI 2013 were utilized to evaluate the performance of the proposed algorithm. The results demonstrate that the segmentation algorithm achieved the best performance with an SSM of 9500 nodes, which obtained a dice of 0.907 and an average surface distance of 1.85 mm. Compared with other methods, our algorithm delineates the prostate region more accurately and efficiently. In addition, we verified the impact of model elasticity augmentation and the fine-tuning item on the network segmentation capability. As a result, both factors have improved the delineation accuracy, with dice increased by 10% and 7% respectively. Conclusions: Our segmentation method has the potential to be an effective and robust approach for prostate segmentation.

as it is crucial for the clinical diagnosis, therapeutic procedure, and treatment planning of various prostate disorders (e.g., prostate cancer, prostatitis or prostatic hypertrophy) b2 2 .

For example, prostate delineation is widely applied for the precise localization of prostate boundary in radiotherapy for dose planning. Besides, in the image-guided computer-assisted surgery, the segmentation of the prostate on preoperative MRI is an essential reference for the inter-operative low-quality image, like ultrasound image b3 3 . However, until now prostate on MR images is still mostly segmented manually by radiologists. The handcrafted delineation of prostate boundary is a time-consuming and labor-intensive operation with a low reproducibility because of its high dependence on medical experience. Moreover, those problems are further aggravated when the borderline is indistinct. As pointed out by Yu et al. b4 4 , automatic prostate segmentation is also a challenging task due to the issue of intensity inhomogeneity, variation of anatomical appearance, and lack of boundary discriminability.

In order to address this challenging task, different automatic or semi-automatic segmentation approaches have been reported in recent years. Martin et al. b5,b6 5,6 proposed a semi-automatic prostate segmentation method, in which a rigid intensity-based registration algorithm and a non-rigid hybrid registration framework were employed successively to align an atlas to the patient image. In their work, 18 MRI series were involved to construct the atlas. Two accuracy metrics respectively based on volume and surface distance were used to investigate the segmentation performance. Results showed that the segmentation accuracy of the apex region and the central region is higher than the base part. In another publication, to add additional knowledge into the segmentation procedure, Korsager et al.

b7 7 combined the spatial information of a prostate atlas with the intensity information in a graph cut segmentation framework to achieve automatic prostate delineation. Their validation experiment was investigated on 76 axial MR images. As a result, a mean Dice similarity coefficient (DSC) of 0.88 and a mean surface distance of 1.45 mm were reported.

Besides, Tian et al. Recently, the performance of deep learning approaches has outperformed to the tranditional state-of-the-art methods in many fields, especially in computer vision b9,b10 9,10 and medical image processing

II. Materials and methods

:section2

Fig.

figure1

1 is an overview of the proposed prostate segmentation algorithm. The segmentation network comprises two branches, serving as a boundary predictor and a label classifier respectively. As shown by Fig. figure1 1A, the first SSM-Net branch is used to deduce from the input, three variables: the SSM global transform, shape control parameters, and the pointwise fine-tuning vector. Then, a patient-specific deformation field is calculated based on the predicted variables, and a surface of target is further generated via a coordinate sampler.

The second ResU-Net branch is employed to predict a probability label map from the input images, instead of directly using its binarized result as the target region. The parameters of the whole network can be optimized by minimizing the loss function value, which is related to two parts (Sec.

sec::section2.4

II.D.): 1) the dice coefficient between the predicted probability map and the input binary label, defined as Loss 1. 2) the similarity between the SSM deformed surface and the input prostate, defined as Loss 2. In the inference phase(Fig. 1st author name or however authors to be briefly identified generated contour obtained by deforming the SSM by using the deformation field predicted by the SSM-Net branch, a distance map is directly calculated to serve as the target boundary constraints. Then, the weighted sum of the distance map and the probability map is regarded as the final prostate segmentation.

N Number of SSM nodes M Number of SSM variations ū3N×1 (ū 3N ) SSM mean shape ψ 3N ×M SSM variation / PCA basis κ M ×1 (κ M ) Variance of SSM variations θ M ×1 (θ M ) SSM shape control parameters t 3 : T, R
Transform parameters: translation and rotation I 3×3 (I) Identity matrix

B. Operators and functions

diag( v) Diagonal matrix of vector v V (ψ 3N ×M , κ M , θ M , ū3N )
Model surface deformation

According to the theory of the PCA-based model, an arbitrary shape can be represented by superimposing a deformation field to the mean shape. In our work, a deformation field is described as the sum of global transform, weighted variations, and the mean shape. That means, according to the notations defined in Table Table1 1, shape can be written as follows:

u = ū3N + V (ψ 3N ×M , κ M , θ M , t 3 ) = ū3N + ψ 3N ×M • diag(κ M ) • θ M • R + T (1) equation1
The deformation ability of SSM depends on its node number (N ) and the variation matrix ψ 3N ×M . As only very limited datasets were used to build the SSM, the model is insufficient to explain all possible shape variations. Two approaches were employed to solve this problem. Firstly, according to the principle of statistic shape, augmenting the example Last edited Date :

II.A. Statistical shape model 1st author name or however authors to be briefly identified shapes by involving small and very smoothly varying deformations, can make the variation matrix (noted as

ψ 3N ×M • diag(κ M ) in Eq. ( equation1 1)) more representative. ψ 3N ×M • diag(κ M ) is
PCA dimensionality reduction form of variation matrix ψ 3N ×M .

ψ 3N ×3N =    k SSM (x 1 , x 1 ) • • • k SSM (x 1 , x N ) . . . . . . . . . k SSM (x N , x 1 ) • • • k SSM (x N , x N )    (2) equation2 k SSM (x, x ) = 1 n -1 n i=1 (u i (x) -µ SSM (x)) (u i (x ) -µ SSM (x )) T (3) equation3
Combining sample covariance kernel (k SSM (x, x )) calculated from sample data, and Gaussian kernel k (s,σ) g can enlarge the flexibility of the model. For the Gaussian kernel in this paper, the smoothness σ = 10mm and the scale of the deformation s = 2mm.

k aug (x, x ) = k SSM (x, x ) + k (2,10) g (x, x ) (4) equation4 k (s,σ) g (x, x ) = s • diag(exp(- x -x 2 σ 2 )) (5) equation5
Secondly, supplementing a point-by-point item ξ 3N to the deformed target (u ) was another feasible approach to represent more possible targets, as shown as Eq. (

6). Item ξ 3N equation6
was predicted by the SSM-net branch automatically. As shown in Figure 2, the prediction of the offset vector shared the same residual CNN structure as the prediction of shape control parameters (noted as θ M in Eq. (

)).

u = ū3N + V (ψ 3N ×M , κ M , θ M , t 3 ) + ξ 3N = ū3N + ψ 3N ×M • diag(κ M ) • θ M • R + T + ξ 3N (6) equation6
The surface model is commonly defined in the physical spatial coordinate system to maintain the visualization invariance on different platforms. The transformation between model space and image space in this work is shown as follows: control parameters with the size of M × 1 (M is the number of shape variations). The involved convolution layers following the "contracting path" with the kernel size 3 × 3 × 3 used 1 pixel stride, and the employed max-pooling layers used pooling size of 2. Similarly, the prediction of fine-tuning vector shared the same residual CNN structure, expect that the average pooling layer is replaced by a max-pooling layer with the pooling size of 2, the stride of 2. In the last step of the SSM-net branch, the variables including the transform matrix, the shape control parameters, and the fine-tuning vector, are input to the last spatial transformation layer, to yield the prostate contour.

P i = (P m -P 0 )/s + 0.5 (
As illustrated by the ResU-net branch in Fig. figure2 2, a residual U-Net is employed to infer the probability label map. U-net b24 24 is a widely used network with high accuracy for object segmentation. In our segmentation framework, each residual block consists of two convolution layers with a kernel size of 3 × 3 × 3, pixel stride of 1. And the max-pooling layers use a pooling size equal to the stride and the size of the up-sample layers is set to 2.

The sigmoid function is utilized as the activation function of the last layer to limit output values to [0,1]. For the whole network, the structure of each layer is shown in Fig.

II.C. Grid mapping ection2.3

As shown in the overview of the segmentation framework (Fig. figure1 1), the patient-specific prostate shape can be obtained by superimposing the predicted deformation field to the standard SSM surface in the inference procedure. While in the training phase, in order to calculate the loss, a predicted binary surface image for SSM (annotated as "generated model surface" in Fig. figure1 1.) is generated by interpolating the input boundary based on the deformation field.

The output binary surface of SSM

g t ∈ R L t W t H t is defined on a regular grid G t = {G t i } = {(x t i , y t i , z t i )}, i ∈ [L t W t H t ]
, where L t ,W t ,H t represent the length, width and height of the output. Similarly, let

G s = {(x s i , y s i , z s i )}, x s i ∈ [0, L s ], y s i ∈ [0, W s ], z s i ∈ [0, H s ]
be the input grid, where, L s ,W s ,H s are the length, width and height of the input binary mask respectively.

The relationship between the output grid G t and the input grid G s can be written as follows.

II. MATERIALS AND METHODS

II.C. Grid mapping figure3

For ∀i ∈ [1...L t W t H t ], G s i = D i (G t i ) =   x t i + d i,x y t i + d i,y z t i + d i,z   =   x s i y s i z s i   (8) equation8 D i (G t i )
is the deformation field predicted by SSM-Net branch, with size of

L t × W t × H t × 3. Based on Eq. ( equation8
8), for each node in the output grid, we can trace its corresponding position on the input mask according to the deformation field. It means that the gray value of output point G t i = (x t i , y t i , z t i ) depends on the gray value of the relevant position

G s i = (x s i + d i,x , y s i + d i,y , z s i + d i,z
) in the input mask.

II.D. Loss function

ection2.4

The loss function of the proposed network consists of two parts: the part for SSM-Net branch and the part for ResU-Net branch.

For SSM-Net branch, according to the mean shape of SSM and the predicted deforma- 1st author name or however authors to be briefly identified tion field, we can calculate the final segmentation directly. To evaluate the accuracy of the prediction, the specific energy function of SSM-Net branch was employed to optimize the deformation field, which is described as follows:

E = L(g t G s ) + λ ϕ = L(g t G s ) + λ 1 θ + λ 2 t + λ 3 ∇ξ (9) equation9
where g represents the gray label in image space. Specifically, g t denotes the boundary of the input mask, and (θ, t, ξ) is the predicted deformation field. In this paper, we defined L as:

L(g t G s ) = 1 -(g t G s )/N (10) equation10
Then, the optimization function can be written as follows:

ϕ = (θ, t, ξ) = arg min(1 -(g s G t )/N ) + λ 1 θ + λ 2 t + λ 3 ∇ξ (11) equation11
where N represents the node number of SSM.

To obtain g t G s , an interpolation is required to calculate the gray value of arbitrary position G s i = (x s i , y s i , z s i ). In this work, 3D bilinear interpolation method b25,b26 25,26 is employed.

Let g t i represents its gray level at point G t i = (x t i , y t i , z t i ). g t i G s i mean the gray level at point G s i = (x s i , y s i , z s i ), and {(m, n, p)} s represents the point (m, n, p) of input mask respectively. According to Eq. (

4), For ∀i ∈ [1...L t W t H t ], g t i G s i = g s (x s i ,y s i ,z s i ) = L s m W s n H s p g t (m,n,p) • max(0, 1 -|x s i -n|) • max(0, 1 -|y s i -m|) • max(0, 1 -|z s i -p|) (12) equation12 equation4
The partial derivatives with respect to gray g t (m,n,p) and coordinate position (x t i , y t i , z t i ) for the backpropagation of loss can be written as follows (∂g

t i G s i /∂y s i /, ∂g t i G s i /∂z s i are similar with ∂g t i G s i /∂x s i ): ∂g t i G s i ∂g t (m,n,p) = L s m W s n H s p max(0, 1 -|x s i -n|)• max(0, 1 -|y s i -m|) • max(0, 1 -|z s i -p|) (13) equation13
II. MATERIALS AND METHODS II.D. Loss function

∂g t i G s i ∂x s i = L s m W s n H s p g t (m,n,p) • max(0, 1 -|y s i -m|) • max(0, 1 -|z s i -p|) •        0, |x s i -n| > 1, -1, x s i > n, 1, x s i < n. ( 14
) equation14
For ResU-Net branch, the following loss function is adopted, where, S mask represent the ground truth. ϑ is the network parameters and P ResU is the output probability map. The values of range from 0 to 1. In the prostate inference step, the binarized weighted sum of the deduced probability map and distance map is regarded as the final prostate segmentation. In terms of probability label map, a bigger value means a higher probability for a pixel to belong to the prostate region. The distance map obtained from the boundary image predicted by the SSM-branch, is such that only the pixel on the boundary has a value of 1, while others are equal to 0. It is calculated according to the following equation: 1st author name or however authors to be briefly identified distance between P i and Pi . 10 is the calculation range, which should be changed according to the image resolution because it decides the prostate voxel range in images. As the image volume was resampled to the same resolution, the calculation range is a constant in our work. Furtherly, as the interior prostate gland being segmented should be evaluated with a large value to reflect its high probability, for the pixel P i in the interior region with D ssm (P i ) less than 0.5, its D ssm value is reassigned to 1.

L(ϑ, S mask ) = 1 - 2 × P ResU (ϑ) × S mask P ResU + S mask (
D ssm (P i ) = 1 -P i -Pi /10 (
The probability map generated by the Res-Unet branch is such that the central region of the prostate has a probability of 1 whilst the marginal prostate region is valued with a lower value, which contributes to most of the prediction deviation, especially when the contour of the prostate is indistinct. In contrast, in this case, the SSM-Net branch can deduce a relatively reasonable boundary due to the representation of prostate shape prior and provided complementary information for the result of the ResU-Net branch. However, the distance map calculated from the SSM-Net branch might not match the segmentation results from the ResU-Net branch, therefore, we have investigated the segmentation performance of the proposed method under different combinations of SSM-Net and ResU-Net branches, as shown as follows:

P combined = w 1 P ResU + (1 -w 1 )D ssm (19) equation19
Two metrics were introduced to evaluate the performance of the proposed segmentation framework, including the dice similarity coefficient (DSC) and the average over the shortest distance between the boundary points of the volumes(ABD). The DSC is formulated as follows:

Dice = 2 * S pred ∩ S mask /( S pred + S mask ) (20) equation20
where S pred and S mask respectively present the predicted segmentation and the input ground truth.

III. Results

:section3

III.A. Data acquisition and experiment set up 1st author name or however authors to be briefly identified

). The experiments showed that the best result in terms of the Dice coefficient is achieved when λ 1 , λ 2 , and λ 3 are equal to 0, 0.01, and 0.01. In the following work, λ 1 , λ 2 and λ 3 are set to 0, 0.01, and 0.01. In addition, we figured out the proposed segmentation framework performs best when the "hyperparameter" w 1 in Eq. ( 5 illustrates the DSC and ABD of the proposed prostate segmentation framework when adopted different SSMs. 4-fold cross-validation was conducted for each group. The DSC and ABD of the ResU-Net branch are constant over the different groups. In terms of the SSM-Net branch, when the network adopts SSMs with 1625 to 9750 nodes (referred as network nodes in the following text), the dice result has significant improved from 0.69 to 0.90 (paired t-test, p <0.001). At the same time, the ABD value has an opposite steep trend, decreasing from 2.63 mm to 2.39 mm. Both dice and ABD have the best result when network node number is 9750, with an average dice of 0.862 and an ABD of 2.04 mm. As regards the performance of the whole framework, it has a similar trend with the SSM-Net branch, with the dice reached a peak of 0.907 and the ABD declined to the lowest of 1.85 mm. For the two groups with more than 9750 nodes, they performed slightly inferior with a Last edited Date : III.C. Accuracy evaluation and analysis page 16

1st author name or however authors to be briefly identified figure5 dice of 0.89. According to the record data, we concluded that the SSM with 9750 nodes is optimal to employ in our framework for the representation of the prostate spatial boundary.

The segmentation procedure of the proposed framework is shown in Fig. figure6 6. The purple model in (c) represents the SSM, whose center is initially positioned at the origin (0,0,0) of the anatomical coordinate system of LPS (Left, Posterior, Superior). As shown in (c), according to the input image, the SSM-Net branch respectively predicts the global transform, weight parameters, and an offset vector for the calculation of the deformation field. The 3D white surface in (d) represents the generated prostate boundary by applying the deformation field to the SSM. Fig. the exhibition, we deduced that the segmentation on the prostate central zone has higher accuracy than the base of the prostate (more complex to delineate). For the case of the severe hyperplastic prostate gland which is larger than the mean shape of SSM, the segmentation framework has relatively poor performance with the maximum ABD was 2.7 mm and the dice coefficient was 0.83. Compared to the segmentation approach employing only the SSM method, our framework can achieve more accurate delineation when the target is beyond the SSM deformation range.

Last edited Date :

III.D. The influence of network flexibility on segmentation accuracy According to the analysis in Section 1st author name or however authors to be briefly identified As shown in Figure 4, driven by their first three principal components of variations, the first four statistical shape models with nodes number from 1625 to 9750 perform quite different from each other, while the two statistic shape models with nodes number of 13000 and 15250 have almost similar performance with the model with 9750. Thus, the mesh shape with nodes number of 9750 is accurate enough to represent the anatomical structure of a regular prostate in the physical coordinate system. However, when nodes number is more than 9750, the accuracy (including DSC and ABD) does not improve with the increase of nodes number. That's because when the model owns a substantial number of nodes, one pixel in the image space may correspond to more than one node in the physical space, resulting in decreased performance.

In the past three decades, three major categories of automatic prostate segmentation The proposed method has a very good potential for clinical application. After establishing the SSM model and training the network, the MR images containing the prostate can be segmented automatically. The 3D model of the prostate can be reconstructed from the segmented results, which is useful for morphological analysis, volume calculation, etc. Moreover, such workflow can be integrated into our previously developed image-guided surgical system b39,b40 39,40 to improve its efficiency and automation.

V. Conclusions

:section5

In this study, we introduced a novel registration-based algorithm that combines CNN and SSM and applied it to the task of precise prostate segmentation. A two branches structure was designed, through which the prior knowledge introduced by SSM and boundary features

V. CONCLUSIONS

extracted by the CNN were fully used for prostate segmentation. Extensive experimental results conducted on two public datasets demonstrated that the proposed network can achieve better performance than several state-of-the-art algorithms for prostate segmentation.

  Materials and methods II.A. Statistical shape model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4 II.B. CNN architecture of the segmentation framework . . . . . . . . . . . . . . . 7 II.C. Grid mapping . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 II.D. Loss function . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 11 12 II.E. Inference of prostate region . . . . . . . . . . . . . . . . . . . . . . . . . . . III. Results III.A. Data acquisition and experiment set up . . . . . . . . . . . . . . . . . . . . . 12 III.B. SSM establishment and analysis . . . . . . . . . . . . . . . . . . . . . . . . . 13 III.C. Accuracy evaluation and analysis . . . . . . . . . . . . . . . . . . . . . . . . 14 III.D. The influence of network flexibility on segmentation accuracy . . . . . . . . Comparison with other methods . . . . . . . . . . . . . . . . . . . . . . . . . (MR) imaging becoming an increasingly important non-invasive imaging modality b1 1 , prostate MR segmentation has been paid close attention in recent years,

  b8 8 utilized a superpixel-based graph cut framework to acquire the prostate surface on MRI. A superpixel is a group of pixels that have similar characteristics such as intensity or location. Due to the capacity to carry wider information, serving as a more Last edited Date : convenient and compact representation of the original image, the superpixel image has been widely used in image segmentation algorithms. In Tian's study, a graph cuts algorithm and an active contour model were integrated for cross-promotion. According to their experiment results, the verification on 43 MRI examples obtained a mean dice of 0.893.

  figure1

1B), according to the Last edited Date : page 4

 4

2 Figure 2 :

 22 Figure 2: The structure of the whole network.

figure2 2 .

 2 figure2 2. The right column of

  figure2

Figure 3 :

 3 Figure 3: The schematic diagram of segmentation generation and image interpolation. (a) shows that the spatial structure of a predicted segmentation is equal to the sum of the SSM and a deformation field. (b) illustrates the 2D calculation strategy of the gray value on the model boundary according to the binary input mask and the deformation field. indicated with red arrows. g(P) represents the gray value of point P.

  15) equation15Thus, the parameterized ResU-Net branch can be optimized during the training procedure: the loss function of the whole network is defined as Eq. ( equation17 17), to optimize the network weights. Where α and β are trainable parameters. L combined = αL(l, M, φ) + βL(ϑ, S mask ) (17) equation17 II.E. Inference of prostate region ection2.5

ection3. 1

 1 We validated the proposed network on two public datasets: the MICCAI PROMISE12 challenge dataset b27 27 and NCI-ISBI 2013 challenge dataset b28 28 . The first database contains 50 III. RESULTS prostate transversal T2 MRIs for training and 30 prostate images for testing. And the second database respectively involves 60, 10, and 10 cases for training, leaderboard, and testing. The two datasets share 11 common volumes, and the ground truths of the PROMISE12 testing set are unavailable. Thus, a total of 119 T2 MRI image cases are collected for our experiment. As the gathered data have different voxel spacing and image size, the input images were resized to the shape of 176 × 128 × 48 beforehand. From the 119 volumes, 40 randomly selected cases contributed to the SSM establishment, and the rest data were used for the network optimization. Specifically, 63 out of the 79 image volumes (3/4) were applied for the network training and the remaining (1/4) for the validation. The modeling procedure and the network segmentation performance are illustrated in Section sec::section3.2 III.B. and Section sec::section3.3 III.C. respectively. And Section sec::section3.4 III.D. analyzes and discusses the influence of SSM flexibility on the framework segmentation accuracy. The network training and all experiments were conducted on a computer with Intel® Core (TM) I7-8700K with a 3.70 GHz CPU, 8 GB memory, and two graphic cards of 8GB NVIDIA GeForce GTX 1080 Ti.III.B. SSM establishment and analysisection3.2The node number and the principal component number are two dominant factors for the flexibility of SSM. The former is decided according to the node number of the counted meshes, and the latter is determined by the compactness of SSM. The compactness is measured via the accumulation of SSM variations which are arranged according to their eigenvalues. In this work, we selected the first M principle components to keep 95% of the total eigenvalues.To investigate the influence of the two factors on the segmentation accuracy of the segmentation method, six SSMs with various nodes and principal components were built. 3D slicer (https://www.slicer.org) was utilized for the refinement of 3D surfaces. In addition, the Gaussian process model building, model fitting, and the PCA-based model building are completed via the Statismo library b22 22 . The details of SSMs are presented in Table

Table2

  2. Thethird column shows the models after flexibility augmentation. Last edited Date : III.B. SSM establishment and analysis page 14

Figure figure4 4

 4 Figurefigure44 illustrates the flexibilities of the primary SSMs. The green models in the middle column depicts the mean shapes with increasing node numbers from top to button. From the middle to the rightmost, the three columns respectively exhibit the deformed shapes generated from the mean shape of SSM by weighting the first three principles variations with three times of their corresponding deviations. Accordingly, the three-column shapes on the left describe the generated models deformed by negative triple deviations. Generally, each row horizontally reveals the influence of the first three principal variations on the deformation of each SSM, and each column vertically compares the different performances of SSMs. In conclusion, the last three SSMs with node numbers of 9750, 13000, 15250 had similar interior and exterior characteristics, while the first three SSMs behaved quite differently.

equation19 19 )

 19 equals 0.6 by trial and error. III. RESULTS III.C. Accuracy evaluation and analysis

Figure 4 :

 4 Figure 4: The flexibility of the statistical shape models with different node numbers. The middle column with green color exhibits the mean models, and the left three columns and the right three columns respectively show the deformed models drove by 3 √ λ times of the first three principal components of variations. λ is the corresponding deviation of each component.

  figure4

Fig.

  Fig.

  figure5

Figure 5 :

 5 Figure 5: The DSC and ABD results of the whole framework, SSM-Net branch, and ResU-Net branch respectively.

figure77

  figure7 7 illustrates the delineation results of four images series. Four slices range from number 16 to number 46 with an interval of 10, are selected to display the recognition

Figure 6 :

 6 Figure 6: The whole segmentation procedure. (a) The input MR image for segmentation. (b) The probability map. (c) The purple model represents the SSM, whose center initially positioned at the origin (0,0,0). (d) The generated contour of the prostate from the SSM. (e) The distance map. (f) The segmentation result.

Figure 7 :

 7 Figure 7: The segmentation results of four image cases via the proposed method. Each row stands one example, and each column shows their segmentation results of the same layers of different examples.

sec::section2. 2

 2 figure8 8 presents the DSC and ABD values of the SSM-Net branch with various nodes and different utilization situations of augmentation and offset item. As shown by those statistical trend lines, regardless of the model augmentation or fine-tuning item, the SSM-Net branch with 9750 nodes or more outperformed other situations. In terms of the model flexibility augmentation, its application has improved the

Figure 8 :

 8 Figure 8: The influence of the deformation ability of SSM-Net branch on its segmentation performances, including dice scores and ABD values). The network elastic ability is mainly determined by three factors: node numbers, flexibility augmentation and the offset item.

figure8

  figure8

  figure8 8, compared with the results of the network without augmentation and offset, the best average dice and distance of the network adopting model elastic augmentation were improved to 0.81 and 2.36 mm respectively. And the network employing both augmentation and offset had higher average dice and smaller average distance than the network utilizing only offset item, as informed by the lines with green and red dots in the figure. Similarly, we figure out that the employment of offset item also contributed to the improvement of the SSM-Net branch. Specifically, compared with the segmentation result of the network without augmentation (the line with cyan dots in subgraph (a)), the involvement of fine-tuning item (the line with blue dots in subgraph (a)) increased the highest average dice of the six groups from 0.74 to 0.79. Based on the network which applies only flexibility augmentation (the lines with blue dots), the utilization of offset items (the lines with red dots) improved the average dice coefficient by 0.06 and the distance by 0.25 mm, making the SSM-Net branch achieve the best result.

  Several teams have combined deep learning algorithms with atlas or deformable models in their researches. Similarly, our segmentation combined neural network with SSM for high precise MRI prostate segmentation. In this way, prior knowledge is introduced by SSM to serve as the boundary constraint and a rough reference, and furtherly combined with a neural network (ResU-Net in our work) which can obtain details from the target image, for the precise recognization of the target region. Besides, as SSM is built based on lots of medical image data, the generated model can reasonably represent the deformation. Therefore, our method required less time (approximately 3 s) while performed a satisfactory segmentation accuracy with high robustness.In medical image processing, the collection of the training dataset limits the application of various learning-based approaches. Fortunately, SSM utilizes geometric information rather than intensity information. Therefore, all images in different modalities such as computed tomography, ultrasound and MRI can contribute to the construction of SSM. It is worth mentioning that, if the dataset is sufficient, active shape model(ASM) and statistical deformation models (SDM) could become superior training supervisors than SSM, as the former can carry intensity information of images and the latter is capable to represent the statistical information of deformation field of a collection of examples. In addition, the introduction of the finite element model (FEM) is worth considering for our further development, because the biomechanical information in FEM can contribute to the delineation of a specific target.

  b11,b12 11,12 . In those researches, neural networks often work as information extractors to eliminate the tedious procedure of traditional feature choice and collection. For precise segmentation of MRI prostate, Guo et al.

	Statistical shape model (SSM) is a geometric model containing a mean shape and mul-
	tiple compressed primarily shape variations of a collection of similar shapes. Due to the
	ability to represent prior geometric information, SSM has been widely applied in different
	medical modalities for the segmentation and registration of various anatomical structures,
	including brain, bone, liver, heart, prostate and so on	b18 18 -	b21 21 . In this paper, based on a
	registration approach, we proposed a novel segmentation algorithm to tackle the prostate
	extraction problem by combining a boundary predictor and a label classifier. Specifically,
	a GoogLeNet-based branch (SSM-Net) was involved as the boundary predictor to deduce
	the prostate contour to obtain a boundary distance map, serving as the target contour con-
	straints. Meanwhile, a 3D residual U-net branch (ResU-Net) was employed as the label
	classifier to predict a probability label map from the input images, to judge the class pos-
	sibility of each pixel. In the inference step, the optimal weighted sum of the distance map
	and the probability map was regarded as the final prostate segmentation. In our validation
	experiment, six different SSMs with various nodes were built to investigate the accuracy and
	efficiency of the whole algorithm.

b13

13 used more concise and effective hierarchical features from MRI prostate image by utilizing a stacked sparse autoencoder. Based on the extracted features, a sparse patch matching method was employed to deduce the corresponding prostate likelihood map, which was further combined with a sparse shape model for the final segmentation. Besides, Mun et al. b14 14 integrated encoding, bridge, decoding, and classification modules to develop a baseline convolutional neural network (CNN) to extract volumetric information. In the meantime, Jia et al. b15 15 researched a coarse-to-fine algorithm for MRI prostate segmentation through a deep learning method. In their algorithm, a registration-based segmentation was firstly used to obtain a rough prostate region. Then a pixel-wise recognizer based on a neural network was further adopted to classify the prostate boundary from the image background. Finally, a refinement algorithm was applied to smooth the contour. Similar to Jia's work, He et al. b16 16 exploited another coarse-to-fine prostate segmentation system via different algorithms. They firstly proposed an adaptive feature learning probability boosting tree for prostate pre-segmentation. Next, a CNN method was developed to obtain the prostate profile model by the judgment of the inner, external, and boundary points. For the last step, an active shape model was employed for the final surface refinement. Their results showed the method is accurate and robust for prostate segmentation, as the neural network was utilized for the extraction of latent image features and the prediction of the prostate boundary. In addition, Wang et al. b17 17 introduced a 3D deep-supervised full CNN with group dilated convolution, aiming to preserve extra image information for prostate delineation. Their method achieved a dice of 0.86. Generally, the aforementioned methods demonstrated that compared with the traditional segmentation algorithms, the approaches based on the deep learning method can delineate the target more accurately in less time. However, the aforementioned ways only involved the information of I. INTRODUCTION the specific input image while without any prior knowledge constraints that are potentially helpful to improve the prostate segmentation accuracy.

Table 1 :

 1 Mathematic notations of statistical shape model

		Figure 1: The overview of the proposed prostate segmentation frame-
	figure1	work.
		II.A. Statistical shape model
	ection2.1
		Due to the capacity to carry prior geometric information of numerous examples in differ-
		ent medical modalities, SSM has been widely applied in object recognition, image process,
		surgery implant design	b22 22 . Generally, it involves two parts to describe the statistical spatial
		information of a collection of objects: a geometric model for the representation of the mean
		shape, and a series of variation vectors to depict the principal components of divergences
		between the objects and the mean shape. As the most prevalent SSM type, Principal Com-
		ponent Analysis (PCA) based shape model can model the variability of various types of
		objects such as images, displacement fields, surface meshes, and volumetric meshes. In our
		II. MATERIALS AND METHODS	II.A. Statistical shape model

Table1

A. Variables: scalars, vectors, matrices

Table 2 :

 2 The details of statistical shape models with various number of nodes

	Table2
	Number of nodes	Original SSMs	Augmented SSMs
		Variation(3N x M) Variation(3N x M)
	1625	4875 × 49	4875 × 50
	3250	9750 × 49	9750 × 50
	6500	19500 × 49	19500 × 50
	9750	29250 × 49	29250 × 50
	13000	39000 × 48	39000 × 50
	16250	48750 × 47	48750 × 50
	III.C. Accuracy evaluation and analysis
	ection3.3

Table 3 :

 3 Comparison between our method and other automatic segmentation methods. All the methods were trained and tested on the same dataset. appearances, and indiscernible boundaries. Therefore, in this work, we proposed a novel prostate segmentation framework, based on CNN and SSM, which has been widely used in prostate segmentation b29 29 . The results demonstrate that the network has the optimal dice of 0.907 and DSC of 1.85 mm under the network nodes of 9750. And both model elastic augmentation and offset applications have positive effects. The performance on the collected clinical data demonstrates that our prostate segmentation framework is feasible, and it has the potential to be a useful clinical tool for the diagnosis, treatment design, and therapeutic procedure of variable prostate disease.

	Table3
	Work	Method	DSC + std ABD [mm]) Time
	Maan et al. Vincent et al. b32 32 b33 33 Toth et al. b36 36 Ou et al. b30 30 Gao et al. b31 31 Milletari et al. b35 35 V-Net 3D AAM AAM Deformation landmark AAM Multi-atlas Multi-atlas + patch-based voxel 0.82 ± 0.02 2.86 ± 0.82 30 min 0.81 ± 0.13 3.09 ± 0.96 4 min 0.86 ± 0.07 2.17 ± 0.63 8 min 0.77 ± 0.18 3.64 ± 1.39 3 min 0.84 ± 0.06 2.85 ± 0.72 40 min 0.86 ± 0.11 2.13 ± 0.86 <1 min Yu et al. b4 4 Volumetric ConvNet 0.87 ± 0.24 2.05 ± 0.69 <1 min Karimi et al. b34 34 CNN + SSM 0.88 ± 0.09 2.16 ± 0.77 <1 min Tian et al. b37 37 PSNet 0.86 ± 0.40 2.72 ± 0.90 <1 min Jia et al. b38 38 Multi-atlas +VGG-19 0.92 ± 0.05 1.63 ± 0.38 40 min
	Ours	SSM-Net + ResU-Net	0.90 ± 0.08 1.85 ± 0.75 <1min
	III.E. Comparison with other methods
	ection3.5
	We compared the performance of our method with existing automatic segmentation methods
	(including atlas-based method, deformable model-based method, and deep learning-based
	method). Non-open-source methods were excluded for comparison. All deep-learning-based
	methods are trained on the same dataset (consisting of 63 volumes), and tested on the
	same dataset (consisting of 16 volumes). For methods that are not based on deep learning,
	their performances are evaluated on the same test dataset (consisting of 16 volumes).. The
	comparison results are summarizes in Table	Table3 3. Results illustrated that our method can
	achieve the second best accuracy than other methods. The inference time consumption of
	our method is much less than atlas or deformable model based methods, and is comparable
	with deep learning network methods.
	IV. Discussion
	:section4
	Prostate segmentation facilitates the diagnosis and treatment of prostate diseases. For ex-

ample, the determination and location of the prostate region are essential information for radiotherapy or high-intensity focused ultrasound operations. However, its clinical application is still limited, because of the segmentation challenges like inhomogeneous intensity, IV. DISCUSSION III.E. Comparison with other methods various anatomical

  have been introduced, including atlas-based algorithms, deformable model-based approaches, and deep learning-based methods. For the first category, the atlas-based strategy has been widely utilized in medical image segmentation and registration b30 30 -b31 31 . Its major principle is to align an atlas that contained spatial prior knowledge to the target images by registration approaches. Then apply the alignment information to deform the atlas label for the final target segmentation. Secondly, in terms of the model-based method, a deformable model is firstly constructed for the representation of standard prostate contour. Then the information extracted from the target image was further applied to drive the model to generate the specific shape. Several groups in the list b32 32 -b34 34 employed deformable model (AAM, ASM, SSM) for the prediction of prostate. For the third group, deep learning-based approaches, especially CNNs are widely introduced in medical image processing because of the powerful feature extraction and non-linear learning capability b34,b4,b35 4,34,35 .

Last edited Date :

Acknowledgments

This work was supported by grants from the National Natural Science Foundation of China (81971709; M-0019; 82011530141), PHC CAI YUANPEI Program(41366SA), Agence Nationale de la Recherche (Investissements d'Avenir Program -grants ANR-11-LABX-0004 and ANR-19-P3IA-0003), Région Auvergne Rhône-Alpes (ProNavIA project), the Foundation of Science and Technology Commission of Shanghai Municipality (19510712200; 20490740700), Shanghai Jiao Tong University Foundation on Medical and Technological Joint Science Research (YG2019ZDA06; YG2021ZD21; YG2021QN72; YG2022QN056), and Hospital Funded Clinical Research(21XJMR02), Xinhua Hospital Affiliated to Shanghai Jiao Tong University, School of Medicine.

Conflict of Interest

The authors have no relevant conflicts of interest to disclose.