An Invasive Mussel ( Arcuatula senhousia , Benson 1842) Interacts with Resident Biota in Controlling Benthic Ecosystem Functioning

: The invasive mussel Arcuatula senhousia has successfully colonized shallow soft sediments worldwide. This ﬁlter feeding mussel modiﬁes sedimentary habitats while forming dense populations and e ﬃ ciently contributes to nutrient cycling. In the present study, the density of A. senhousia was manipulated in intact sediment cores taken within an intertidal Zostera noltei seagrass meadow in Arcachon Bay (French Atlantic coast), where the species currently occurs at levels corresponding to an early invasion stage. It aimed at testing the e ﬀ ects of a future invasion on (1) bioturbation (bioirrigation and sediment mixing) as well as on (2) total benthic solute ﬂuxes across the sediment–water interface. Results showed that increasing densities of A. senhousia clearly enhanced phosphate and ammonium e ﬄ uxes, but conversely did not signiﬁcantly a ﬀ ect community bioturbation rates, highlighting the ability of A. senhousia to control nutrient cycling through strong excretion rates with potential important consequences for nutrient cycling and benthic–pelagic coupling at a broader scale. However, it appears that the variability in the di ﬀ erent measured solute ﬂuxes were underpinned by di ﬀ erent interactions between the manipulated density of A. senhousia and several faunal and / or environmental drivers, therefore underlining the complexity of anticipating the e ﬀ ects of an invasion process on ecosystem functioning within a realistic context.


Introduction
Globalization, translated into increasing the transport of goods and people, aqua-and agriculture practices, contributes to increasing the risk of species introductions to new areas across the world [1]. Estuarine and coastal ecosystems concentrate on activities associated with the above-mentioned vectors of introduction and are also facing high levels of several anthropogenic stressors (e.g., eutrophication, habitat destruction, chemical pollution). Thus, they constitute hot spots for biological invasions [2]. Invasive species, while exhibiting a rapid spatial and temporal expansion of their population, may profoundly modify biodiversity-ecosystem functioning relationships (BEF) in their host ecosystem [3,4]. This can clearly affect ecosystem services delivered to human well-being in estuarine and coastal areas [5]. Assessing the effect of introduced species appears of particular interest for conservation and management. One main issue consists of the sometimes considerable time lag between the first introduction of the species to the spread of its population [6] (i.e., before the population actually becomes invasive). As the species during the lag phase exists in very low numbers in habitats where it has never been recorded, it might be undetected for years, overlooked, or even misidentified [7]. Studying the potential effects of the species on the surrounding environment and on ecosystem function at this early stage may, however, help target future management actions [8,9].
Numerous publications on biological invasions have reported the effects on the structure of the invaded community [10][11][12][13]. Some studies have found a decline in native species populations due to competition or predation [14,15]. Others have found facilitative effects through habitat provision in the form of physical structures [16,17]. Some studies have also focused on their effects on ecosystem functioning [3,4,[18][19][20][21][22]. Such effects are dependent on both the traits of the invasive species [23] and properties of the ecosystem they invade [16,24] and are, therefore, by default, context-dependent [25]. Due to this context-dependency, ecological effects can also change in direction and magnitude following changes in environmental conditions, for example, in [26], or being enhanced when occurring in systems under other intense anthropogenic pressures [27]. Therefore, understanding the impact of the invader and the relative importance of other potential drivers on the populations of the invader itself and native species and on the studied ecosystem functions are of importance for a proper assessment of the ecological consequences of the invasion.
The Asian date mussel or bag mussel, Arcuatula senhousia (Benson, 1842) (Mytilidae), is a native of East Asia from Siberian coasts to Singapore. It has been accidentally introduced to the Pacific coast of North America, Australia, New Zealand, the Mediterranean Sea, and the French Atlantic coast [7,[28][29][30]. The occurrence of A. senhousia has also recently been recorded in UK waters [31], the Netherlands [32], West Africa [33], and in the Azov-Black Sea Basin [34]. The mussel inhabits intertidal and shallow subtidal soft sediments [7,28,30,35]. A short life cycle, rapid growth, and a prolonged reproductive period result in large population fluctuations [28,36,37], but these opportunistic traits also contribute to enhancing their establishment and spread within colonized areas [28,38]. As for other filter-feeding bivalves such as oysters whose development into reefs actually drive ecosystem functioning and organic matter cycling in estuaries and coastal bays [39][40][41], A. senhousia may form dense mats where individuals attach together through (the formation of) byssus threads, thereby strongly altering sedimentary habitats [42]. The development of such mats indeed seems to be conditioned by mussel density, environmental conditions, and/or the local biodiversity of the resident fauna [43]. Densities up to ca. 15,000 individuals per square meter have been reported in southern California or in New Zealand [44], and up to 5000 individuals in Mediterranean lagoons [38].
Various effects of A. senhousia have been reported from both the native and invaded areas. In some areas, they serve as prey for many organisms such as fishes, birds, crabs, or predatory gastropods [45][46][47]. Through their important excretion rates, they efficiently contribute to nutrient cycling [42,48], sediment organic enrichment [49], and are thus considered as carbon dioxide generators, potentially increasing the transfer of CO 2 to the atmosphere [38]. Magnitude and direction of the effects of the organic enrichment and biodeposition seem to be partly density-dependent, mostly depending upon the development of byssal mats than can be detrimental by smothering the seagrasses that the species often co-occurs with [50] and destroying local fisheries [43]. It, however, also creates a 3D structure habitat promoting abundances of some species through facilitation [42,51]. Interestingly, A. senhousia was found associated with the red algae Agarophytion vermiculophyllum in its native range before major shifts in ecosystem structure, leading to extended decline of seagrass meadows [43], likely having a habitat cascading effect together with other ecosystem engineers such as seagrasses and/or red algae [52]. In low densities, they can facilitate deposit-feeding organisms [30,51]. A study found evidence of the mussels feeding on a common food source than the other suspension feeding bivalves Cerastoderma glaucum, Ruditapes decussatus, and Scrobicularia plana [53], hence possibly promoting interspecific competition for resources [54]. Competitive aptitudes of A. senhousia could then lead to the limitation of carbon transfers from phytoplankton to native consumers [20]. Such an effect of A. senhousia is clearly further enhanced by a rapid recolonization of that particular species in soft sediment following dramatic disturbance events such as the development of anoxic conditions [36].
In Arcachon Bay (France), A. senhousia was first reported in 2002 and has since been found in various locations in the Bay (Figure 1), being present in Zostera noltei meadows, bare sediments, and oyster (Magallena gigas) reefs [7,55,56]. This is clearly indicative of a successful colonization at the scale of the Bay. Up to now, A. senhousia in Arcachon Bay has, however, never reached densities comparable to those recorded in highly invaded ecosystems nor formed dense continuous mats, with a maximum density of 22 indiv.m −2 recorded in 2007 [7]. It has nevertheless been pointed out that there was a significant contribution of A. senhousia to community biomass in sparse Z. Noltei meadows where macrofauna density was low [55] and recent observations in April 2019 reported slightly higher adult A. senhousia density of 100 indiv·m −2 [57] within a comparable habitat. Colonization of the Bay by A. senhousia then appeared to still be limited, although this situation may correspond to an early stage of invasion process preceding a boom in population corresponding to an exponential growth [6]. Therefore, studying the effects of A. senhousia on biological community and ecosystem functioning at this stage is of particular interest in order (1) to assess the current effect of this non-native species within its host ecosystem as well as (2) to anticipate potential consequences of a future population boom, and in the context of other anthropogenic disturbances affecting benthic habitats such as the multifactorial strong decline of Z. noltei meadows observed during the last 30 years [58]. Mussel population establishment could indeed affect the activity and/or the persistence of other organisms in the habitat they invade through various density dependent physical, biological, or biogeochemical processes [30,42,46,51].

General Strategy
Experiments were performed from 27 March to 11 April 2019 at the marine station of Arcachon. It consisted in the addition of an increasing number of A. senhousia individuals in natural sediment cores containing resident macrofauna community and seagrasses. These cores were used to perform lab-based bioturbation and nutrient flux measurements through incubations in a tidal system.

Animal and Sediment Cores Collection
A total of 227 A. senhousia individuals were manually sampled at the study site on 27 March 2019 and brought back to the lab. Their maximum shell length was measured using calipers and placed in a flow-through system until being transferred to the experimental cores.
A total of 27 PVC tubes (h = 30 cm, internal diameter = 9.5 cm) were pushed down in the sediment and retrieved, and the height of the sediment interface was adjusted so it just matches with the top of the tube using foam stoppers from underneath. These were then brought back to the lab and placed in the tidal mesocosm system (see below). In the present study, we aimed at looking at the consequences on benthic ecosystem functioning of a potential future invasion of A. senhousia leading to an increase in population density, within a sparse seagrass meadow where the species currently occurs at levels likely corresponding to an early invasion stage. Hence, by manipulating ex situ A. senhousia densities in intact sediment cores, we allowed for the assessment of mussel density effects and their interactions with the resident biota while measuring proxies of macrofauna community activities (bioturbation including sediment particle mixing and bioirrigation) and of ecosystem functioning (solute fluxes). We aimed to disentangle the relative importance of the density-dependent effects of A. senhousia on bioturbation and solute fluxes across the sediment-water interface from those of the resident biota during an invasion process and their potential consequences.

Study Site
Arcachon Bay is a meso-to macrotidal lagoon located in the French Atlantic coast (Figure 1). It is directly connected to the Bay of Biscay and mostly receives freshwater inputs from the Leyre River located in the eastern side. The major part of intertidal flats is colonized by the dwarf eelgrass Zostera noltei, thereby constituting what is considered as the largest European intertidal seagrass meadow [59]. The experiment was based on samples taken at the "Cassy" site (44 • 42 24.4" N, 1 • 03 38.4" W) in the eastern end of Arcachon Bay (SW France, Figure 1). The site was characterized by the presence of a sparse Z. noltei meadow occupying a muddy intertidal flat and characterized by a mosaic of low shoot-density seagrasses and bare sediment patches [60].

General Strategy
Experiments were performed from 27 March to 11 April 2019 at the marine station of Arcachon. It consisted in the addition of an increasing number of A. senhousia individuals in natural sediment cores containing resident macrofauna community and seagrasses. These cores were used to perform lab-based bioturbation and nutrient flux measurements through incubations in a tidal system.

Animal and Sediment Cores Collection
A total of 227 A. senhousia individuals were manually sampled at the study site on 27 March 2019 and brought back to the lab. Their maximum shell length was measured using calipers and placed in a flow-through system until being transferred to the experimental cores.
A total of 27 PVC tubes (h = 30 cm, internal diameter = 9.5 cm) were pushed down in the sediment and retrieved, and the height of the sediment interface was adjusted so it just matches with the top of the tube using foam stoppers from underneath. These were then brought back to the lab and placed in the tidal mesocosm system (see below).

Experimental Setup
The experimental setup consisted of one main 300 L tank where the cores were disposed upright. Tidal cycles were reproduced using a second 300 L tank that was used to transfer the water every six hours, resulting in a 12 h long tidal cycle. Light was provided (114 µmol photon m −2 ·s −1 ) from above the cores during low tides occurring during the day. Cores otherwise received ambient outdoor light. During the course of the experiment, temperature and salinity of the water used for the incubation varied between 11.9 (night time) and 15.1 • C (day time), and between 32.1 and 33.0 • C, respectively.

Experimental Treatments
The following day, 0, 4, 8, 12 or 24 A. senhousia individuals were added in the cores, corresponding to the control, low, medium, high, and very high density treatments, respectively (five replicates per treatment). Mean shell length per core was comprised between 24.5 ± 2.6 mm and 27.8 ± 2.0 mm and did not significantly differ amongst cores (one-way univariate PERMANOVA). During experimental low tide conditions, mussels were gently pushed into the sediment to ca. 1/3 of their total length. Most of the mussels quickly buried and only shew their siphons out of the sediment while other individuals tended to stay at the sediment surface and/or aggregated ( Figure 2). They were thereafter left for four days acclimatization time before the start of the incubations.

Bioturbation Measurements
Bioirrigation and particle mixing were respectively quantified in each core through the measurement of the depth profiles of solute (bromide) and particle (luminophores) tracers after 10 days of incubation. Br-concentration was set at 10 mM in the tank using a NaBr solution before the start of the incubation and 2.5 g dry weight (DW) of green luminophore tracers (Environmental Tracing EcoTrace ® , Helensburgh, UK, D 50 = 35 µm, density = 2.5 g·cm −3 ) were suspended, homogenized in seawater, and gently spread at the sediment core surface using a Pasteur pipette at the same time. Water was constantly aerated using air-stones during the whole course of the incubation. Temperature and salinity were monitored using a Hobo sensor placed at the surface of a dedicated core.
At the end of the experiment and for each core, a small sample of sediment was taken by gently scraping the core surface for carbon and nitrogen content analyses. Porewater was subsequently sampled using rhizon pore water samplers inserted in holes fitted on the side of each core at 0.5, 1.0, 1.5, 2.0, 3.0, 4.0, 5.0, 6.0, and 7.0 cm depth in the sediment. Br-concentrations in each sample were then determined spectrometrically [61] in order to generate vertical Br-profiles thereafter fitted with the enhanced-diffusion model. In this case, solute transport induced by the biological activity is quantified with a factor of enhanced diffusion (ε) [62]. Cores were then finally sliced at 0.5, 1, 1.5, 2, 3, 4, 5, 6,7,8,9,11,13, and 15 cm depths. Slices were homogenized and an approximatively 30 g aliquot of sediment sampled for luminophore counting after ensuring that no macrofauna were trapped. The rest of the core was sieved through a 1 mm sieve for further macrofauna sorting and identification. Sediment aliquots were freeze-dried and 1 g of dry sediment photographed under UV light using a digital camera. Luminophore pixels were counted after a binarization step (based on the RGB level) for each image using image analysis software [63]. Vertical depth profiles were then computed based on the relative concentrations of luminophores in each slice and fitted using a continuous time random walk (CTRW) model [64,65]. So-obtained normal biodiffusion coefficient (D b N in cm 2 ·yr −1 ) values reflect the particle mixing intensity by resident macrofauna [18,66].

Nutrient Flux Measurements
One day prior the end of the experiment, nutrients (Si, NO x , NH4 + , PO4 3− ) and dissolved oxygen total fluxes between the sediment and the water column were measured in each core. Cores were retrieved from the tank and fitted with 0.5 L incubation chambers filled with experimental tank water and provided with magnetic stirrers. They were thereby incubated in darkness for two hours. At the beginning and the end of the incubation, water samples were taken in each core for Si, NO x , NH 4 + , and PO 4 3− concentration measurements while oxygen total fluxes were obtained through the linear fitting of the decrease in oxygen saturation measured using oxygen micro-optodes (firesting minioptodes, OXROB10, Pyroscience ® , Aachen, Germany) 2 points calibrated using air saturation as 100% and a sodium ascorbate solution set to 0%.

Macrofauna and Seagrass Characteristics
Sieved material (>1 mm) obtained from the cores was fixed in 4% formalin and macrofauna (including A. senhousia individuals) was then sorted, identified under stereomicroscope at the lowest taxonomic level, counted and weighed after being dried at 60 • C for 48 h.
Above-and below-ground biomasses of seagrasses were obtained in each core after being separated and dried at 60 • C for 48 h.

Statistical Analyses
The aims of the analyses were to: (1) assess the effects of the tested increasing density of A. senhousia on community bioturbation rates as well as benthic oxygen and nutrient total fluxes, and (2) contextualize such effects with regard to the characteristics of the invaded ecosystem. We then compared the amount of variations imputable to A. senhousia to those directly linked with the variability in habitat characteristics and the other species/functional groups naturally present in the community.
Differences in main environmental characteristics (sediment C and N contents and ratios, Z. noltei above-and below-ground biomass) among density treatments were tested through one-way univariate permutation analysis of variances (PERMANOVAs) (factor TR, 5 levels) based on Euclidean distance and using 999 permutations. Pairwise PERMANOVAs were performed as well as associated dispersion analyses (PERMDISP) in order to highlight differences in data variability.
The contribution of A. senhousia compared to the other members of the macrofauna community as well as environmental factors in explaining the variability in the measured bioturbation metrics, and in solute fluxes were investigated using a distance-based redundancy analysis (dbRDA) performed with the DistLM option in the PERMANOVA + add-on for PRIMER [67]. Forward selection was used to build models using AIC selection criterion. To reduce the number of predictor variables and to eliminate species with only one entry in the data that could complicate further analysis, we grouped the species according to their functional attributes. Analyses were carried out on the basis of biomasses. This resulted in seven macrofauna predictors (treatment-wise densities and dry weights are found in Table 1): A. senhousia dry weight, dry weight of burrowing polychaetes (including Glycera tridaclyla, Hediste diversicolor, and Nepthys hombergii), S. plana dry weight, dry weight of filter feeding bivalves (Including Cerastoderma edule and Ruditapes philippinarum), Abra segmentum dry weight, and Tubificoides benedii dry weight. Peringia ulvae was excluded from the statistical analyses since they were free to move between the cores and their presence in the cores at the time of incubation could thus not be verified. Environmental predictors included Z. noltei above-and below-ground DW as well as surface sediment C and N contents. To equalize the multivariate variances for further analyses, the dry weight of burrowing polychaetes was square root transformed and the dry weights of filtering bivalves and S. plana were log(x + 1) transformed. Table 1. Variables measured in experimental cores (treatment-wise mean, maximum, minimum and standard deviation, n = 5): macrofauna abundance (AB) and dry weight (DW) (Arc: A. senhousia, Filt biv: filtering bivalves, Burr poly: burrowing polychaetes, Sp: S. plana, As: A. segmentum, Tb: T. benedii, Pu: P. ulvae), environmental variables (C/N ratio of surface sediment, dry weight of above-and belowground Z. noltei in g·m −2 ), bioturbation (Normal biodiffusion coefficient D b N for particle mixing and biodiffusion enhancement factor ε for bioirrigation) and fluxes across the SWI.

Results
The recovery percentage for A. senhousia was 98.3 % overall and ranged from 95 to 100 % between the treatments. Resulting tested densities are presented in Figure 3A and

Bioturbation
The Db N significantly differed between treatments (Permanovas F' = 4.0629, p = 0.018). Such result was more likely to be related with high intra-treatment variability, and thus especially in the medium treatment for which a maximum Db N individual (core) value of 20.08 cm²·yr −1 was measured ( Figure  4, Table 1), which led to the detection of a significant difference compared to the control. The

Environmental Variables
Experimental cores were randomly sampled at the collection site within a large Z. noltei meadow patch, only avoiding invasive red algae Agarophyton vermiculophyllum (formerly Gracilaria vermiculophylla) that was apparently forming a complex with the seagrass and A. senhousia by a possible facilitation mechanism [52]. Such sampling procedure led to inter-core variations in environmental variables (Z. noltei above-and below-ground biomasses, sediment C and N content, Table 1). There were, however, no significant differences among treatments in the amount of neither above-nor below-ground biomass of the seagrass. The C and N contents of the surface sediment did not differ between treatments. The mean C and N contents were highest in the low-density treatment and lowest in the medium-density treatment, and tended to increase toward the very high-density treatment but the differences were not statistically significant.

Solute Fluxes
Phosphate and ammonium fluxes across the water-sediment interface were significantly affected by A. senhousia density (main PERMANOVAs tests gave F' = 9.5171 with associated p = 0.002 and F' = 3.1674 with associated p = 0.01 for phosphate and ammonium fluxes, respectively). Fluxes were lowest in the control treatment (phosphate 0.34 ± 0.93, ammonium 0 ± 0 µmol·m −2 ·h −1 , mean ± SD throughout) and increased with density showing significant differences compared to the control treatment from low density treatment for phosphate fluxes and from medium density treatment for ammonium fluxes ( Figure 3C,D). Mean fluxes reached their maximum values in the very high density treatment (phosphate 13.29 ± 4.95, ammonium 137.5 ± 149.7 µmol·m −2 ·h −1 ) ( Figure 3; Table 1). Fluxes of phosphate and ammonium moreover correlated strongly (r = 0.79).
Although no global effect of A. senhousia density on oxygen fluxes across the sediment-water interface was detected (PERMANOVA, F' = 2.1187, p = 0.12), they were, however, significantly higher at the high and very high densities compared to the control ( Figure 3B) to reach values of −3.34 and −4.16 mmol·m −2 ·h −1 (Table 1, Figure 3B).
NOx and Si fluxes were not affected of A. senhousia and were also characterized by high intra-treatment variability (Table 1, Figure 3E,F).

Bioturbation
The D b N significantly differed between treatments (Permanovas F' = 4.0629, p = 0.018). Such result was more likely to be related with high intra-treatment variability, and thus especially in the medium treatment for which a maximum D b N individual (core) value of 20.08 cm 2 ·yr −1 was measured ( Figure 4, Table 1), which led to the detection of a significant difference compared to the control. The enhancement biodiffusion coefficient did not significantly differ between treatments (Permanovas F' = 0.5754 with p = 0.689) as well as with a well-marked intra-treatment variability ( Figure 4, Table 1). enhancement biodiffusion coefficient did not significantly differ between treatments (Permanovas F' = 0.5754 with p = 0.689) as well as with a well-marked intra-treatment variability ( Figure 4, Table 1).

Contextualization of A. senhousia Effects
To further investigate the relative importance of faunal and environmental drivers behind the tested effect of A. senhousia density, we conducted separate variation partitioning procedures, DistLM (PRIMER PERMANOVA+) with each solute flux and measured the bioturbation metric as the

Contextualization of A. senhousia Effects
To further investigate the relative importance of faunal and environmental drivers behind the tested effect of A. senhousia density, we conducted separate variation partitioning procedures, DistLM (PRIMER PERMANOVA+) with each solute flux and measured the bioturbation metric as the response variables and macrofauna as well as environmental variables (sediment C and N content, above-and below-ground biomass of Z. noltei) ( Table 1) as the predictor variables.
Not surprisingly, the biomass of A. Senhousia was the most important driver for phosphate and ammonium fluxes, respectively accounting for 58 and 37% of their variability. A total of 76% of the variance of phosphate fluxes could moreover be significantly explained when including the biomasses of S. plana and A. segmentum in the model in addition to the biomass of the mussels (Table 2), while a model including the above-ground biomass of seagrasses, those of S. plana and of burrowing polychaetes in addition to the biomass of mussels accounted for 55% of the variance in ammonium fluxes ( Table 2). Table 2. Sequential distance-based linear models result between predictors (macrofauna and environment) and response variables (benthic fluxes across the sediment water interface (PO 4 3− , NH 4 + , NO x , Si, and O 2 ) and community bioturbation rates (particle mixing normal biodiffusion coefficient D b N and bioirrigation biodiffusion enhancement factor ε)). Models were identified using the forward selection procedure and AIC selection criterion. Bold font indicates significant results (p < 0.05). Oxygen fluxes, according to the DistLM results, conversely seemed to be primarily driven by the aboveground (i.e., leaves) of seagrasses. However, A. senhousia accounted for 18% of the oxygen flux variance within a model also including the biomass of S. plana, all together explaining 61% of the total oxygen flux variations.

Response
Arcuatula senhousia biomass only accounted for a small amount (9%) of the variability in Si fluxes, mainly driven by S. plana. Mussel biomass did not account for any variation in NO x or the bioturbation metrics. Instead, the D b N pattern was clearly driven by S. plana and A. segmentum (12 and 7% of variation accounted for, respectively), though these results were not statistically significant due to high variability in the occurrence of the species in the replicate cores. The variation in biodiffusion enhancement coefficient (ε) was significantly accounted for by burrowing polychaetes and Z. noltei below-ground biomass, and these two variables positively correlated with the ε coefficient. According to the results from the DistLM analyses, A. senhousia accounted for 58 and 37% of the variation in the phosphate and ammonium fluxes, respectively, and was clearly the most important predictor of these fluxes. For phosphate, other large-sized bivalves, S. plana and A. segmentum, also significantly directly affected the flux (12 and 6%, respectively).

Discussion
Non-native species invasion processes are controlled by many intertwined biotic and abiotic factors. These include the physiological and ecological characteristics of the species as well as the natural inherent properties (physical, biological) of the invaded ecosystem [16,23,50] that are also potentially altered by anthropogenic pressures [27]. Therefore, forecasting the level at which a given invasive ecosystem engineer would actually drive the ecosystem functioning at the scale of the host habitat (while overwhelming the effects of the resident biological community) requires assessing the invader density-effects within a realistic framework (i.e., that encapsulates the variations in habitat structure naturally occurring in the host ecosystem). The density effect of A. senhousia on bioturbation and solute fluxes was assessed during the present study in natural sediment cores where the abundance of the resident fauna taxa and the density of seagrasses varied regardless of the tested density treatments (Table 1). Inter-core variability thereby reflected small scale spatial heterogeneity typically observed within the studied sparse Zostera meadow [60,68,69]. Use of such an experimental design generally led to high intra-treatment variability, somehow challenging to handle from a statistical point of view [18] but clearly answered to an identified need for the incorporation of environmental complexity in view of better assessing the effects of non-native ecosystem engineers on BEF relationships [70]. It also aimed at contextualizing some of the observed effects of the targeted invasive ecosystem engineer assuming that species-specific expression of functional traits could be modulated by changes in the abiotic and biotic environment [18,71,72].
In the present study, we showed that increasing densities of A. senhousia within a sparse Z. noltei meadow clearly enhanced solute fluxes across the sediment-water interface (SWI), especially through increasing phosphate and ammonium effluxes, but conversely did not significantly affect community bioturbation rates. These results are in good accordance with earlier studies highlighting that, besides its fast metabolism and high excretion rates [38,48], A. senhousia has a very limited effect on deep dwelling macrofauna species such as large bivalves [73] known to profoundly control bioturbation rates [74], especially within Z. noltei meadows [66]. Thus, we logically identified these large burrowing bivalves (including S. plana and the filtering bivalves, e.g., the common cockle C. edule and the manilla clam R. philippinarum) and large burrowing polychaetes as the main drivers of the variability in bioturbation rates (Table 2). Large bivalves were also found as the main drivers for the silicate fluxes across the SWI, conversely to A. senhousia only showing a marginal and non-significant effect ( Table 2). Silicate fluxes across the SWI are known to be controlled by bioturbation [75] and biogenic silica dissolution, the later particularly enhanced in bivalve biodeposits containing diatom frustules [21,76]. Thereby, it appears fully coherent that the same faunal drivers (i.e., the main bioturbating invertebrate taxa) were identified as explaining variability in both silicate fluxes and bioturbation rates (sediment mixing and bioirrigation). It therefore suggests that the strong increases in phosphate, ammonium fluxes, and to a lesser extent, in oxygen uptake, measured amongst density treatments are likely related with the fast metabolism and strong excretion rates of A. senhousia, estimated as amongst the largest in mytilids [48], rather than from a triggered sediment organic matter remineralization by the microbial compartment [77]. This is further supported by previous seasonal estimations of biogenic fluxes of phosphate and ammonium across the SWI induced by A. senhousia excretion in tidal flats of its native area (Japan) [78], where the values calculated in spring were in the same order of magnitude than our actual flux measurements at medium, high, and very high density treatments. Corresponding field density of A. senhousia reported in spring by these authors (ca. 1000 to 1500 indiv·m −2 ) being directly comparable with our tested densities. At the relatively short-time scale of our experiment, changes in ecosystem functioning induced by the Asian date mussels then seemed to be linked with increasing community metabolism rather than changes in sediment-fauna interactions. These changes, moreover, could not result from the exclusion of some taxa or plant disappearance as a response to the establishment of A. senhousia and its byssal network [30,50].
Enhancement of community metabolism induced by increasing density of A. senhousia was here revealed through increasing phosphate, ammonium, and oxygen fluxes ( Figure 3, Table 2). Although these density effects on the different solutes supporting ecosystem functioning were quite similar and consistent with each other, it appears that they were underpinned by different interactions between the manipulated density of A. senhousia and several faunal and/or environmental drivers, therefore underlining the complexity of realistic context assessments of the effects of invasive species on BEF relationships.
Significant increase in phosphate fluxes across the SWI were already observed at low density treatment compared to the control (Figure 3c) to reach a maximum value of 13.29 ± 4.95 µmol·m −2 ·h −1 at very high density ( Figure 3, Table 2). These phosphate fluxes were not affected by the presence of seagrasses in the cores, leading to an amount of 76% of variance that could be explained when adding the biomasses of S. plana and A. segmentum to the manipulated A. senhousia biomass ( Table 2). It suggests that establishment of A. senhousia is likely to enhance the recycling of dissolved phosphorus in the benthic compartment and its release in the water column, enhancing benthic-pelagic coupling. Such a boosted phosphorus efflux could be of major importance since phosphate ions are directly used by micro and macrophytes for their growth [79,80] and would thereby support an enhanced primary production by either microphytobenthos at the sediment surface or phytoplankton in the water column. Conversely, resulting impacts on phanerogams such as Z. noltei are more difficult to apprehend. On one hand, Z. noltei leaves are able to uptake phosphorus in the form of dissolved phosphate [81], suggesting that increasing phosphate input to the water column could support their growth. On the other hand, it has been shown at the same study site than ours, that phosphorus is indeed preferentially taken up by the roots system of Zostera from the porewater [60]. Accordingly, early spring new root development is likely to promote oxygen penetration into the sediment, thereby leading to storage of phosphorus through adsorption on iron hydroxides and resulting in a large influx of dissolved phosphate from the water to the sediment during the productive season [60], the so-constituted stock of phosphorus being used as a reservoir for plant growth. Our increasing phosphorus fluxes oriented toward the water column induced by A. senhousia density, even at low density, could therefore be indicative of a change in the phosphorus cycle and the accumulation of phosphorus in the sediment rooted zone, and this added to the fact that A. senhousia strongly inhibited Zostera rhizome growth [50]. This hypothesis is further supported by the increasing oxygen consumption induced by the mussels, likely limiting its penetration in the sediment column and associated phosphorus binding with iron hydroxides. Available phosphorus for seagrass growth that could be taken up by the roots would likely diminish through time, potentially affecting plant development and the seasonal dynamics of nutrient recycling.
Ammonium fluxes measured during the present study were negligible in the control and the low density treatments for the three remaining density treatments ( Figure 3, Table 1). Low ammonium concentrations in the water column have frequently been reported in Arcachon Bay [82][83][84] as well as lower NH4 + fluxes across the SWI measured in Zostera meadows compared to unvegetated areas that were assumed to be the result of inorganic nitrogen uptake by plant roots, rhizomes, and leaves counteracting NH4 + effluxes induced by animal excretion and microbial-mediated organic matter remineralization [85]. Here again, the difference in ammonium fluxes between treatments indeed underlines the fast metabolism and nutrient recycling capacity of A. senhousia compared to those of the other main macrofauna species such as the manilla clam R. phillipinarum [48,78], even at temperatures at which the remaining resident fauna were not able to induce significant fluxes at the scale of our incubations, or to overwhelm the NH4 + uptake by the seagrasses. It should also be underlined that the levels of ammonium flux we measured at sufficient density to do so (i.e., from 846 to 3400 indv·m −2 ) were in the same order of magnitude as the NH4 + fluxes measured at the same study site within the seagrass meadow through dark benthic chamber incubations in May 2006 [60] (i.e., before A. senhousia was present), and at much higher water temperatures (from 20.2 to 21 • C) than in our experiments (ca. 15 • C). Therefore, at the scale of the studied sparse seagrass meadow, an increase in ammonium efflux through mussel excretion could clearly indicate a change in the seasonal dynamics of nitrogen processing with potential consequences, as above-mentioned for phosphate, on the growth or development of micro-or macrophytes at both the sediment surface and the water column including the seagrasses themselves. Indeed, enhanced ammonium effluxes by A. senhousia can participate in supporting Zostera leaf growth [50], but also the development of opportunistic macroalgae that are more efficient at taking up nitrogen from the water column [86]. Accordingly, we observed at the study site during sampling an association of A. senhousia with the red algae A. vermiculophyllum entangled in mussel byssal threads [57]. This highly invasive and widespread macroalgae species is known to be particularly efficient at taking up and thereby removing NH4 + from the water column [86]. It could then be logically hypothesized that A. vermiculophyllum would further benefit from enhanced ammonium fluxes toward the water column induced by A. senhousia metabolism, highlighting probable habitat cascading facilitation processes between the two species, as already observed with other mytilids [52]. Moreover, it has been already demonstrated that A. vermiculophyllum in the Aveiro Lagoon (Portugal) tends to outcompete with Z. noltei, resulting in the long run in a switch from a seagrass meadow to macroalgae dominated mudflats [86]. In Z. Noltei meadows of Arcachon Bay, sediment denitrification processes that could also lead to the removal of nitrogen from the system have been assessed to occur at very low rates regardless of the season [87], and were not measured during the present study. Conversely, it is now established that denitrification rates in mudflats are enhanced by both the establishment of filtering bivalve populations at high densities such as observed during the restoration of oyster reefs [39,40] as well as by the accumulation of A. vermicullophyllum at the sediment surface [88]. It could be hypothesized then that a potential future widespread invasion by both the Asian date mussel and the A. vermiculophyllum macroalgae associated with seagrass disappearance would further modify nitrogen cycling and loss mechanisms for N inputs. It therefore points out the need for further investigations aiming at assessing more precisely the relationships between these two species and their potential impact on the seagrass meadows of Arcachon Bay.
Non-surprisingly, significant inter-treatment differences in oxygen fluxes across the SWI have been detected ( Figure 3) indicating an increase in community respiration rates with increasing A. senhousia density. However, a closer examination at the core scale revealed that changes in above-ground biomass of Zostera amongst cores actually explained a greater part of this variance before the manipulated variance of A. senhousia. In other words, the effect of the Asian date mussel on oxygen fluxes cannot be decoupled in this particular habitat from the Z. noltei leaves biomass and associated respiration rates. The A. senhousia invasion, although enhancing community respiration in the sparse seagrass meadow, then seem to not overwhelm the effects of the respiration by dense (with relatively high biomass) Zostera patches, underlying the context dependency of the effects of invasive ecosystem engineers on ecosystem functioning [18]. Such a complex relationship between the Asian date mussel and seagrasses is in good agreement with the finding of [50], who showed that (1) A. senhousia could affect the growth of Zostera, but also that success of implantation of A. senhousia depended on seagrass density. Modification of community respiration rates in Zostera meadow in Arcachon Bay associated with an invasion of A. senhousia would be of particular concern because of the crucial importance of this biogenic habitat, known to (at least partly) control the primary production and global oxygen dynamics at the scale of the whole Arcachon Bay ecosystem [59,69].
In a scenario of a future expansion of the invasion of sparse Z. noltei meadows by A. senhousia at levels comparable to what has been observed in several ecosystems, it can be expected that ecosystem functioning would be strongly enhanced because of the important excretion rates and fast nutrient recycling capacity of the mussels. It then appears likely that seasonal dynamics of both phosphorus and nitrogen cycling would be likely controlled by the mussel's excretion rates such as observed in its native area [48,78], in contrast to a current nutrient seasonal cycle triggered by seagrass metabolism and microbial-mediated sediment processes [85,87]. Moreover, such changes in nutrient cycling would further facilitate the establishment of a habitat cascade implying seagrasses, mussels, and the macroalgae A. vermiculophyllum, and are likely to affect the persistence of the seagrass meadows [89], already strongly threatened at the scale of the whole Bay [58]. These findings clearly highlight the importance of monitoring the spreading of A. senhousia within Z. noltei meadows, and this already at an early stage of invasion, in order to better assess and anticipate changes in ecological functioning of this key habitat. More precisely, seasonal variations in the spatial distribution of the Asian date mussel is characterized by successive extension and retreat phases in its native area [90]. It then appears of particular importance to assess these spatial and temporal variations over a larger scale within the extended meadows to target future management strategies.