# Increasing densities of an invasive polychaete enhance bioturbation with variable effects on solute fluxes

Laura Kauppi, Guillaume Bernard, Ralf Bastrop, Alf Norkko, Joanna Norkko

► To cite this version:

Laura Kauppi, Guillaume Bernard, Ralf Bastrop, Alf Norkko, Joanna Norkko. Increasing densities of an invasive polychaete enhance bioturbation with variable effects on solute fluxes. Scientific Reports, 2018, 8 (1), 10.1038/s41598-018-25989-2. hal-04568242

HAL Id: hal-04568242

https://hal.science/hal-04568242

Submitted on 4 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
SCIENTIFIC REPORTS

<figure><img src="image_1.png" /></figure>

OPEN

# Increasing densities of an invasive polychaete enhance bioturbation with variable effects on solute fluxes

Received: 20 October 2017

Accepted: 25 April 2018

Published online: 16 May 2018

L. Kauppi ID<fnref n="1" />,², G. Bernard<fnref n="1" />,³, R. Bastrop⁴, A. Norkko<fnref n="1" />,⁵ & J. Norkko ID<fnref n="1" />

Bioturbation is a key process affecting nutrient cycling in soft sediments. The invasive polychaete genus *Marenzelleria* spp. has established successfully throughout the Baltic Sea increasing species and functional diversity with possible density-dependent effects on bioturbation and associated solute fluxes. We tested the effects of increasing density of *M. arctia*, *M. viridis* and *M. neglecta* on bioturbation and solute fluxes in a laboratory experiment. Benthic communities in intact sediment cores were manipulated by adding increasing numbers of *Marenzelleria* spp. The results showed that *Marenzelleria* spp. in general enhanced all bioturbation metrics, but the effects on solute fluxes varied depending on the solute, on the density and species identity of *Marenzelleria*, and on the species and functional composition of the surrounding community. *M. viridis* and *M. neglecta* were more important in predicting variation in phosphate and silicate fluxes, whereas *M. arctia* had a larger effect on nitrogen cycling. The complex direct and indirect pathways indicate the importance of considering the whole community and not just species in isolation in the experimental studies. Including these interactions provides a way forward regarding our understanding of the complex ecosystem effects of invasive species.

Benthic invertebrates have a pivotal role in modifying biogeochemical properties and nutrient transformation processes at the sediment-water interface<fnref n="1" />⁻³. Nutrient remineralization is a key ecosystem process providing the basis for primary production in pelagic and benthic food webs⁴⁻⁶. Remineralization processes are affected by a number of factors including temperature, organic matter quantity and quality, and the structure of benthic infauna (both macro- and micro-organisms), modified by seasonality in the system⁷⁻¹⁰. Macrofauna enhances remineralization processes through their bioturbation, which enhances the activity of the microbial communities in the sediment ultimately responsible for organic matter remineralization²,⁹,¹¹, through oxygenating the sediment which enhances aerobic mineralization¹²,¹³, and by enhancing diffusion<fnref n="1" />. Bioturbation activities by benthic macrofauna include all transport processes carried out by animals that directly or indirectly affect sediment matrices. These processes include both particle reworking and burrow ventilation¹⁴. Differences in species richness and the species' functional traits, as well as interactions among these determines their contribution to ecosystem processes, such as nutrient cycling/remineralization¹⁵,¹⁶ and may be even more important for ecosystem functioning than species richness alone⁷,¹⁷⁻²⁰.

Motivated originally by the effect of species loss on ecosystems, biodiversity-ecosystem functioning (BEF) studies have concluded that increasing biodiversity, and especially functional diversity, enhances ecosystem functioning in experimental conditions²¹,²². Species invasions are feared for their potential negative consequences on biodiversity, even though this fear is often not evidence-based²³. On the contrary, species invasions may even increase the species and functional diversity in some areas²⁴,²⁵, but the importance of this potential increase for ecosystem functioning remains largely unknown²² with further uncertainty about extrapolating effects from short-term experiments to longer-term effects in the ecosystem. Experimental studies on BEF-relationships offer

<footnote n="1">¹University of Helsinki, Tvärminne Zoological Station, J.A. Palménin tie 260, FI-10900, Hanko, Finland. ²Department of Environmental Sciences, P.O. Box 65 (Viikinkaari 1), 00014 University of Helsinki, Helsinki, Finland. ³CNRS, EPOC, UMR 5805, F33400, Talence, France. ⁴University of Rostock, Institute of Biological Sciences, Albert-Einstein-Str. 3, D-18059, Rostock, Germany. ⁵Baltic Sea Centre, Stockholm University, Stockholm, Sweden. Correspondence and requests for materials should be addressed to L.K. (email: laura.kauppi@helsinki.fi)</footnote>
SCIENTIFIC REPORTS

<figure><img src="image_2.png" /></figure>

OPEN

# Increasing densities of an invasive polychaete enhance bioturbation with variable effects on solute fluxes

L. Kauppi ID<fnref n="1" />,², G. Bernard<fnref n="1" />,³, R. Bastrop⁴, A. Norkko<fnref n="1" />,⁵ & J. Norkko ID<fnref n="1" />

Received: 20 October 2017

Accepted: 25 April 2018

Published online: 16 May 2018

Bioturbation is a key process affecting nutrient cycling in soft sediments. The invasive polychaete genus *Marenzelleria* spp. has established successfully throughout the Baltic Sea increasing species and functional diversity with possible density-dependent effects on bioturbation and associated solute fluxes. We tested the effects of increasing density of *M. arctia*, *M. viridis* and *M. neglecta* on bioturbation and solute fluxes in a laboratory experiment. Benthic communities in intact sediment cores were manipulated by adding increasing numbers of *Marenzelleria* spp. The results showed that *Marenzelleria* spp. in general enhanced all bioturbation metrics, but the effects on solute fluxes varied depending on the solute, on the density and species identity of *Marenzelleria*, and on the species and functional composition of the surrounding community. *M. viridis* and *M. neglecta* were more important in predicting variation in phosphate and silicate fluxes, whereas *M. arctia* had a larger effect on nitrogen cycling. The complex direct and indirect pathways indicate the importance of considering the whole community and not just species in isolation in the experimental studies. Including these interactions provides a way forward regarding our understanding of the complex ecosystem effects of invasive species.

Benthic invertebrates have a pivotal role in modifying biogeochemical properties and nutrient transformation processes at the sediment-water interface<fnref n="1" />⁻³. Nutrient remineralization is a key ecosystem process providing the basis for primary production in pelagic and benthic food webs⁴⁻⁶. Remineralization processes are affected by a number of factors including temperature, organic matter quantity and quality, and the structure of benthic infauna (both macro- and micro-organisms), modified by seasonality in the system⁷⁻¹⁰. Macrofauna enhances remineralization processes through their bioturbation, which enhances the activity of the microbial communities in the sediment ultimately responsible for organic matter remineralization²,⁹,¹¹, through oxygenating the sediment which enhances aerobic mineralization¹²,¹³, and by enhancing diffusion<fnref n="1" />. Bioturbation activities by benthic macrofauna include all transport processes carried out by animals that directly or indirectly affect sediment matrices. These processes include both particle reworking and burrow ventilation¹⁴. Differences in species richness and the species' functional traits, as well as interactions among these determines their contribution to ecosystem processes, such as nutrient cycling/remineralization¹⁵,¹⁶ and may be even more important for ecosystem functioning than species richness alone⁷,¹⁷⁻²⁰.

Motivated originally by the effect of species loss on ecosystems, biodiversity-ecosystem functioning (BEF) studies have concluded that increasing biodiversity, and especially functional diversity, enhances ecosystem functioning in experimental conditions²¹,²². Species invasions are feared for their potential negative consequences on biodiversity, even though this fear is often not evidence-based²³. On the contrary, species invasions may even increase the species and functional diversity in some areas²⁴,²⁵, but the importance of this potential increase for ecosystem functioning remains largely unknown²² with further uncertainty about extrapolating effects from short-term experiments to longer-term effects in the ecosystem. Experimental studies on BEF-relationships offer

<footnote n="1">¹University of Helsinki, Tvärminne Zoological Station, J.A. Palménin tie 260, FI-10900, Hanko, Finland. ²Department of Environmental Sciences, P.O. Box 65 (Viikinkaari 1), 00014 University of Helsinki, Helsinki, Finland. ³CNRS, EPOC, UMR 5805, F33400, Talence, France. ⁴University of Rostock, Institute of Biological Sciences, Albert-Einstein-Str. 3, D-18059, Rostock, Germany. ⁵Baltic Sea Centre, Stockholm University, Stockholm, Sweden. Correspondence and requests for materials should be addressed to L.K. (email: laura.kauppi@helsinki.fi)</footnote>
important mechanistic understanding about the ecosystem responses to biodiversity change²⁶, but can also be criticized for the highly controlled conditions used to reduce the natural variability that might confound the results of the mechanistic processes, even if some studies suggest that results from small-scale experiments are generalizable on a larger scale of observational studies²⁷.

Invasive species are capable of a rapid population expansion in both space and time, and they can be either native or non-native. Increasing densities of a single species, however, increases the expression of the functional traits associated with this species, thus possibly affecting the ecosystem functioning depending on the studied process⁵,²⁴. Recent field studies quantifying benthic oxygen and nutrient fluxes in the multiply stressed, highly invaded Baltic Sea have provided evidence of enhanced nutrient cycling in the sediment following the establishment of an invasive polychaete genus *Marenzelleria* spp.¹⁰,²⁸,²⁹. In total three species of the genus have become established and spread into the entire Baltic Sea since the first observations in the southern Baltic in 1985³⁰⁻³², *M. viridis* and *M. neglecta* of North American, and *M. arctia* of Arctic origin³¹. In the study area, *M. arctia* prefers deeper (>20 m), muddy bottoms, whereas *M. neglecta* and *M. viridis* prefer shallower (up to 20 m) bottoms and *M. viridis* especially sandier areas³³. Especially in the deeper, muddy sediments *M. arctia* can be a dominant member of the benthic communities along with the clam *Macoma balthica*, while all three of them have increased species richness and added functionality in the naturally species-poor benthic communities of the Baltic Sea²⁵. Species of the genus *Marenzelleria* are classified as facultative deposit feeders and suspension feeders³⁴. All three species in the Baltic Sea are classified as biodiffusers and have low particle reworking rates³⁵. Compared to the native fauna, *M. viridis* and *M. neglecta* burrow deeper (down to 25–35 cm depth) into the sediment constructing L- or J-shaped burrows, whereas *M. arctia* burrows down to 6–8 cm depth and constructs a network of burrows (ref. ³⁵, pers. obs). *M. viridis* and *M. neglecta* affect the solute transport in the sediment mainly through non-local, advective transport, whereas the solute transport mode of *M. arctia* is of a more diffusive character³⁵. Partly therefore, *M. viridis* and *M. neglecta* are also suggested to affect solute fluxes and nutrient cycling more than *M. arctia*³⁵,³⁶. These conclusions are, however, based on results from highly controlled, single-species experiments and do not necessarily reflect the species' effects in a real-life context with variable community composition possibly affecting the outcome of BEF-relationships.

Phosphorus and nitrogen are the most important macronutrients for primary production, and play a key role in eutrophication. In the long run, increased density of *Marenzelleria* spp. could lead to increased burial of phosphorus in the sediment and thus mitigate eutrophication; increased density of the polychaetes is hypothesized to increase oxygenation of the sediment through bioturbation, which increases the availability of iron-oxyhydroxides that capture phosphate ions³⁷,³⁸. The fate of phosphorus is especially interesting, since the decreased internal phosphorus recycling results in decreased blooms of harmful cyanobacteria, which are phosphorus limited and plague the Baltic Sea ecosystem³⁹. Increased oxygenation of the sediment could also facilitate the return of other macrofauna in previously hypoxic areas. The dominance of the *Marenzelleria* species complex in many areas in the Baltic Sea, and their potential importance for a key ecosystem process, nutrient cycling, also call for a mechanistic understanding of the species' impact in the organic matter remineralization processes.

The effect of macrofauna on solute fluxes depends e.g. on density, biomass and functional traits of individual species and their interactions in the community¹⁷,⁴⁰,⁴¹. In this study, we tested the density-dependence of the effect of *Marenzelleria* spp. on nutrient cycling in a density-manipulation experiment following the rationale in the model of Norkko & Reed *et al*.³⁸. Their results show a release of phosphorus at moderate densities of worms (1000 to 3000 ind m⁻²) and a retention with population densities of over 3000 ind m⁻². Thus, in a simple, single-species model, the direction and magnitude of P flux is dependent of *Marenzelleria* spp. density. In order to have as close a resemblance to natural conditions as possible, we conducted a density-manipulation experiment using intact sediment cores collected from the field. Our aims were to investigate 1) how an increasing density of *Marenzelleria* spp. affects fluxes of nitrate+nitrite (NOₓ) ammonium (NH₄⁺), phosphate (PO₄³⁻) and silicate (Si⁴⁺), 2) what bioturbation metrics are affected by *Marenzelleria* spp. and other macrofauna and 3) which bioturbation metrics are associated with the different solute fluxes?

## Results

**Experimental conditions.** The sediment at the core collection site was classified as sand (D₅₀ = 0.46 mm) with a C/N ratio of 7.05 and C content of 0.007%. The organic matter content of the sediment at the sediment collection site varies between 2.6 and 11.2% (measured as loss on ignition) throughout the year¹⁰. The surface sediment in the cores was a 3–4 cm layer of very fine silt, underneath which there was a 1–2 cm layer of gravel followed again by soft mud. Temperature in the water during the experiment varied between 12.5 to 14.2 °C and salinity varied between 5.5 and 5.8. Oxygen concentration in the cores ranged from 7.44 mg l⁻¹ to 8.62 mg l⁻¹, corresponding to a saturation of 75 to 86% during the experimental period. During the dark incubations the temperature in the incubation tank remained at approximately 14 °C, and salinity was 5.7.

**Marenzelleria spp. and other macrofauna recovered after the experiment.** The community in the control cores represents the background community composition with hydrobid snails and *Marenzelleria* spp. as the most abundant species (see Table 1 for core-wise mean, min and max) followed by *Macoma balthica*, *Cerastoderma glaucum* and *Hediste diversicolor*. The clam *M. balthica* had the largest biomass (measured as g wet weight) in the control cores, followed by *Marenzelleria* spp., *H. diversicolor*, *C. glaucum* and the hydrobid snails. Other species occurring in the experimental cores were Chironomidae, *Pygospio elegans*, *Manayunka aestuarina*, *Potamopyrgus antipodarum*, *Cyanophthalma obscura*, *Oligochaeta*, and some *Mya arenaria*, which could have a very large biomass. Of the *Marenzelleria* spp. added to the density manipulated cores 12.5% to 100% were recovered after the experiment. The low survival (12.5%) in one of the cores was due to large specimens of *H. diversicolor* present and feeding on the *Marenzelleria* spp. Excluding this core the percentage of added worms
important mechanistic understanding about the ecosystem responses to biodiversity change²⁶, but can also be criticized for the highly controlled conditions used to reduce the natural variability that might confound the results of the mechanistic processes, even if some studies suggest that results from small-scale experiments are generalizable on a larger scale of observational studies²⁷.

Invasive species are capable of a rapid population expansion in both space and time, and they can be either native or non-native. Increasing densities of a single species, however, increases the expression of the functional traits associated with this species, thus possibly affecting the ecosystem functioning depending on the studied process⁵,²⁴. Recent field studies quantifying benthic oxygen and nutrient fluxes in the multiply stressed, highly invaded Baltic Sea have provided evidence of enhanced nutrient cycling in the sediment following the establishment of an invasive polychaete genus *Marenzelleria* spp.¹⁰,²⁸,²⁹. In total three species of the genus have become established and spread into the entire Baltic Sea since the first observations in the southern Baltic in 1985³⁰⁻³², *M. viridis* and *M. neglecta* of North American, and *M. arctia* of Arctic origin³¹. In the study area, *M. arctia* prefers deeper (>20 m), muddy bottoms, whereas *M. neglecta* and *M. viridis* prefer shallower (up to 20 m) bottoms and *M. viridis* especially sandier areas³³. Especially in the deeper, muddy sediments *M. arctia* can be a dominant member of the benthic communities along with the clam *Macoma balthica*, while all three of them have increased species richness and added functionality in the naturally species-poor benthic communities of the Baltic Sea²⁵. Species of the genus *Marenzelleria* are classified as facultative deposit feeders and suspension feeders³⁴. All three species in the Baltic Sea are classified as biodiffusers and have low particle reworking rates³⁵. Compared to the native fauna, *M. viridis* and *M. neglecta* burrow deeper (down to 25–35 cm depth) into the sediment constructing L- or J-shaped burrows, whereas *M. arctia* burrows down to 6–8 cm depth and constructs a network of burrows (ref. ³⁵, pers. obs). *M. viridis* and *M. neglecta* affect the solute transport in the sediment mainly through non-local, advective transport, whereas the solute transport mode of *M. arctia* is of a more diffusive character³⁵. Partly therefore, *M. viridis* and *M. neglecta* are also suggested to affect solute fluxes and nutrient cycling more than *M. arctia*³⁵,³⁶. These conclusions are, however, based on results from highly controlled, single-species experiments and do not necessarily reflect the species' effects in a real-life context with variable community composition possibly affecting the outcome of BEF-relationships.

Phosphorus and nitrogen are the most important macronutrients for primary production, and play a key role in eutrophication. In the long run, increased density of *Marenzelleria* spp. could lead to increased burial of phosphorus in the sediment and thus mitigate eutrophication; increased density of the polychaetes is hypothesized to increase oxygenation of the sediment through bioturbation, which increases the availability of iron-oxyhydroxides that capture phosphate ions³⁷,³⁸. The fate of phosphorus is especially interesting, since the decreased internal phosphorus recycling results in decreased blooms of harmful cyanobacteria, which are phosphorus limited and plague the Baltic Sea ecosystem³⁹. Increased oxygenation of the sediment could also facilitate the return of other macrofauna in previously hypoxic areas. The dominance of the *Marenzelleria* species complex in many areas in the Baltic Sea, and their potential importance for a key ecosystem process, nutrient cycling, also call for a mechanistic understanding of the species' impact in the organic matter remineralization processes.

The effect of macrofauna on solute fluxes depends e.g. on density, biomass and functional traits of individual species and their interactions in the community¹⁷,⁴⁰,⁴¹. In this study, we tested the density-dependence of the effect of *Marenzelleria* spp. on nutrient cycling in a density-manipulation experiment following the rationale in the model of Norkko & Reed *et al*.³⁸. Their results show a release of phosphorus at moderate densities of worms (1000 to 3000 ind m⁻²) and a retention with population densities of over 3000 ind m⁻². Thus, in a simple, single-species model, the direction and magnitude of P flux is dependent of *Marenzelleria* spp. density. In order to have as close a resemblance to natural conditions as possible, we conducted a density-manipulation experiment using intact sediment cores collected from the field. Our aims were to investigate 1) how an increasing density of *Marenzelleria* spp. affects fluxes of nitrate+nitrite (NOₓ) ammonium (NH₄⁺), phosphate (PO₄³⁻) and silicate (Si⁴⁺), 2) what bioturbation metrics are affected by *Marenzelleria* spp. and other macrofauna and 3) which bioturbation metrics are associated with the different solute fluxes?

## Results

**Experimental conditions.** The sediment at the core collection site was classified as sand (D₅₀ = 0.46 mm) with a C/N ratio of 7.05 and C content of 0.007%. The organic matter content of the sediment at the sediment collection site varies between 2.6 and 11.2% (measured as loss on ignition) throughout the year¹⁰. The surface sediment in the cores was a 3–4 cm layer of very fine silt, underneath which there was a 1–2 cm layer of gravel followed again by soft mud. Temperature in the water during the experiment varied between 12.5 to 14.2 °C and salinity varied between 5.5 and 5.8. Oxygen concentration in the cores ranged from 7.44 mg l⁻¹ to 8.62 mg l⁻¹, corresponding to a saturation of 75 to 86% during the experimental period. During the dark incubations the temperature in the incubation tank remained at approximately 14 °C, and salinity was 5.7.

**Marenzelleria spp. and other macrofauna recovered after the experiment.** The community in the control cores represents the background community composition with hydrobid snails and *Marenzelleria* spp. as the most abundant species (see Table 1 for core-wise mean, min and max) followed by *Macoma balthica*, *Cerastoderma glaucum* and *Hediste diversicolor*. The clam *M. balthica* had the largest biomass (measured as g wet weight) in the control cores, followed by *Marenzelleria* spp., *H. diversicolor*, *C. glaucum* and the hydrobid snails. Other species occurring in the experimental cores were Chironomidae, *Pygospio elegans*, *Manayunka aestuarina*, *Potamopyrgus antipodarum*, *Cyanophthalma obscura*, *Oligochaeta*, and some *Mya arenaria*, which could have a very large biomass. Of the *Marenzelleria* spp. added to the density manipulated cores 12.5% to 100% were recovered after the experiment. The low survival (12.5%) in one of the cores was due to large specimens of *H. diversicolor* present and feeding on the *Marenzelleria* spp. Excluding this core the percentage of added worms
<table><caption>**Table 1.** Core-wise nutrient fluxes used as response variables, and macrofauna densities and biomasses, and bioturbation parameters used as predictors in the DistLM-analysis. Note that the values for *Marenzelleria* include both the natural community and the worms added as part of the density treatments. Bioturbation parameters $D_b^N$ = biodiffusion coefficient, SR = % of surface reworked, MPD = maximum penetration depth, BI = bioirrigation.</caption><thead><tr><th rowspan="2">Core</th><th colspan="2">NO<sub>x</sub></th><th colspan="2">NH<sub>4</sub><sup>+</sup></th><th colspan="2">PO<sub>4</sub><sup>3-</sup></th><th colspan="2">Si<sup>4+</sup></th><th colspan="2">M. arctica</th><th colspan="2">M. neg+vir</th><th colspan="2">Marenzelleria</th><th colspan="2">C. glaucum</th><th colspan="2">M. balthica</th><th colspan="2">H. diversicolo</th><th colspan="2">Hydrobiidae</th><th colspan="2">D<sub>b</sub><sup>N</sup></th><th colspan="2">SR</th><th colspan="2">MPD</th><th colspan="2">BI</th></tr><tr><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">ind m<sup>-2</sup></th><th colspan="2">ww g m<sup>-2</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th></tr></thead><tbody><tr><td>Ca</td><td>-0.02</td><td>-0.15</td><td>-2.08</td><td>-0.08</td><td>0</td><td>1805</td><td>23.43</td><td>0.00</td><td>155.56</td><td>0.00</td><td>3.65</td><td>0.73</td><td>83.23</td><td>6.50</td><td>14.10</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>Cb</td><td>0.05</td><td>0.07</td><td>0.50</td><td>0.62</td><td>211</td><td>1053</td><td>10.27</td><td>3.27</td><td>0.14</td><td>0.00</td><td>7.55</td><td>0.00</td><td>76.83</td><td>5.50</td><td>32.74</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>Cc</td><td>0.10</td><td>0.38</td><td>3.37</td><td>0.85</td><td>0</td><td>1625</td><td>15.00</td><td>7.38</td><td>24.62</td><td>0.09</td><td>5.13</td><td>5.84</td><td>82.79</td><td>5.50</td><td>41.35</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>Cd</td><td>0.45</td><td>0.46</td><td>86.82</td><td>4.90</td><td>181</td><td>0</td><td>1.93</td><td>0.00</td><td>186.82</td><td>49.69</td><td>5.32</td><td>2.19</td><td>57.56</td><td>8.50</td><td>10.98</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>1a</td><td>-0.01</td><td>0.46</td><td>52.85</td><td>0.29</td><td>812</td><td>812</td><td>21.34</td><td>0.00</td><td>238.50</td><td>3.61</td><td>3.19</td><td>5.84</td><td>68.57</td><td>10.00</td><td>78.89</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>1b</td><td>0.03</td><td>0.25</td><td>18.64</td><td>2.98</td><td>0</td><td>1264</td><td>1.75</td><td>0.00</td><td>6.84</td><td>4.89</td><td>5.05</td><td>1.46</td><td>74.49</td><td>5.50</td><td>36.31</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>1c</td><td>0.07</td><td>0.39</td><td>-4.22</td><td>0.10</td><td>542</td><td>0</td><td>2.71</td><td>9.35</td><td>0.00</td><td>3.77</td><td>5.07</td><td>1.10</td><td>58.86</td><td>8.50</td><td>21.72</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>1d</td><td>-0.02</td><td>0.81</td><td>4.59</td><td>-0.86</td><td>650</td><td>975</td><td>28.03</td><td>27.20</td><td>0.00</td><td>28.23</td><td>5.34</td><td>3.29</td><td>77.01</td><td>6.50</td><td>97.63</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>2a</td><td>0.09</td><td>0.61</td><td>18.89</td><td>0.78</td><td>421</td><td>843</td><td>41.75</td><td>6.61</td><td>0.00</td><td>0.45</td><td>4.68</td><td>0.73</td><td>79.12</td><td>5.50</td><td>95.07</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>2b</td><td>-0.06</td><td>0.14</td><td>8.52</td><td>0.41</td><td>316</td><td>948</td><td>7.17</td><td>2.40</td><td>183.93</td><td>0.00</td><td>2.33</td><td>0.37</td><td>69.50</td><td>5.50</td><td>37.63</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>2c</td><td>0.05</td><td>0.25</td><td>-2.31</td><td>0.04</td><td>1011</td><td>1517</td><td>13.32</td><td>2.33</td><td>0.00</td><td>0.00</td><td>2.69</td><td>1.10</td><td>68.69</td><td>4.50</td><td>22.53</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>2d</td><td>0.01</td><td>0.85</td><td>14.60</td><td>-0.45</td><td>963</td><td>481</td><td>10.79</td><td>2.24</td><td>0.02</td><td>0.00</td><td>3.09</td><td>1.46</td><td>71.63</td><td>12.00</td><td>57.80</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>3a</td><td>0.00</td><td>0.44</td><td>-1.69</td><td>0.17</td><td>1341</td><td>1006</td><td>33.57</td><td>2.17</td><td>0.00</td><td>2.94</td><td>2.94</td><td>2.19</td><td>80.05</td><td>7.50</td><td>14.27</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>3b</td><td>0.02</td><td>0.27</td><td>-1.60</td><td>1.26</td><td>963</td><td>481</td><td>4.95</td><td>1.23</td><td>0.02</td><td>0.00</td><td>1.23</td><td>0.37</td><td>65.73</td><td>3.50</td><td>54.24</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>3c</td><td>-0.02</td><td>0.51</td><td>0.57</td><td>-0.26</td><td>722</td><td>1805</td><td>35.18</td><td>11.99</td><td>124.48</td><td>0.00</td><td>7.91</td><td>1.46</td><td>69.93</td><td>8.50</td><td>45.77</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>4a</td><td>0.08</td><td>0.97</td><td>16.56</td><td>0.32</td><td>2949</td><td>843</td><td>34.01</td><td>8.10</td><td>173.39</td><td>0.00</td><td>5.23</td><td>4.02</td><td>80.06</td><td>5.50</td><td>67.37</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>4b</td><td>0.04</td><td>0.89</td><td>1.74</td><td>1.54</td><td>4198</td><td>1399</td><td>25.72</td><td>3.05</td><td>209.76</td><td>0.00</td><td>2.33</td><td>2.19</td><td>86.64</td><td>7.50</td><td>90.18</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>4c</td><td>0.02</td><td>0.93</td><td>8.11</td><td>-0.07</td><td>2785</td><td>464</td><td>20.63</td><td>3.92</td><td>0.00</td><td>14.10</td><td>0.99</td><td>1.10</td><td>80.98</td><td>10.00</td><td>60.01</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>4d</td><td>-0.02</td><td>1.04</td><td>3.76</td><td>0.93</td><td>2437</td><td>4062</td><td>65.47</td><td>2.92</td><td>0.00</td><td>0.00</td><td>2.40</td><td>4.75</td><td>88.35</td><td>8.50</td><td>131.25</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>5b</td><td>-0.09</td><td>1.39</td><td>12.11</td><td>1.65</td><td>9749</td><td>0</td><td>43.86</td><td>4.39</td><td>14.88</td><td>0.00</td><td>4.04</td><td>1.46</td><td>91.34</td><td>4.50</td><td>140.17</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>5c</td><td>-0.01</td><td>-0.70</td><td>-8.27</td><td>-2.19</td><td>433</td><td>650</td><td>14.58</td><td>42.64</td><td>9.64</td><td>164.87</td><td>3.63</td><td>2.56</td><td>90.30</td><td>8.50</td><td>103.81</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>average</td><td>0.04</td><td>0.49</td><td>11.02</td><td>0.62</td><td>1461</td><td>1049</td><td>21.69</td><td>6.72</td><td>63.27</td><td>12.98</td><td>3.99</td><td>2.11</td><td>74.70</td><td>7.05</td><td>59.71</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>min</td><td>-0.09</td><td>-0.70</td><td>-8.27</td><td>-2.19</td><td>0</td><td>0</td><td>1.75</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.99</td><td>0.00</td><td>25.83</td><td>3.50</td><td>10.98</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>max</td><td>0.45</td><td>1.39</td><td>86.82</td><td>4.90</td><td>9749</td><td>4062</td><td>65.47</td><td>42.64</td><td>238.50</td><td>164.87</td><td>7.91</td><td>5.84</td><td>91.34</td><td>12.00</td><td>140.17</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr></tbody></table>

surviving was between 67 and 100%. The realized densities and biomasses of the most common infauna are presented in Table 1.

**Solute fluxes.** NOₓ fluxes in the experimental cores ranged from -0.092 mmol m⁻² d⁻¹ (high-density core 5b) to 0.446 mmol m⁻² d⁻¹ (unmanipulated core Cd), ammonium fluxes ranged from -0.704 mmol m⁻² d⁻¹ (high-density core5c) to 1.388 (high-density treatment 5b), phosphate fluxes ranged from -0.00827 mmol m⁻² d⁻¹ (high-density core 5c) to 0.0866 mmol m⁻² d⁻¹ measured in one of the controls (Cd), and silicate fluxes ranged from 2.186 mmol m⁻² d⁻¹ (in 5c) to 4.897 mmol m⁻² d⁻¹ (in Cd). Core-wise solute fluxes are presented in Table 1.

**Bioturbation metrics.** At the end of the experiment, burrow networks created by *Marenzelleria* spp. were seen along the core walls. Luminophore tracers were buried down into these burrows (Fig. 1A). Typical burrows created by *H. diversicolor* and siphonal gallery networks generated by *M. balthica* during deposit feeding and associated with tracer burial were also apparent. We could also observe fecal-pellet strings, typical of *Marenzelleria* spp., around the openings of the burrows at the sediment surface (Fig. 1B). These pellets were composed of both sediment and luminophore tracers, indicating that *Marenzelleria* spp. individuals had deposit-fed at the surface during the experiment.

The biodiffusion coefficient $D_b^N$ ranged from 0 (unmanipulated core Cb) to 5.84 cm² d⁻¹ (low-density manipulated core 1a). The percentage of surface reworked (SR) ranged from 57.6% (low-density manipulated core 1c) to 91.3% (high-density manipulated core 5b). The maximum penetration depth (MPD) ranged from 3.5 cm (core 3b) to 12.0 cm (core 2d). Bioirrigation (BI) had a minimum of 10.98 ml d⁻¹ in the unmanipulated core Cd, and a maximum of 140.17 ml d⁻¹ in the high-density treated core 5b. Core-wise bioturbation metrics are presented in Table 1.

**Pore water profiles and ΔNH₄⁺, ΔNOₓ, ΔPO₄³⁻ and ΔSi⁴⁺ during dark incubation.** At the start of the incubation, the pore water profiles of the control core Ca (Fig. 2a,b) showed profiles with increasing concentrations of ammonium, phosphate and silicate, and decreasing concentration of nitrate+nitrite from the surface down. During the incubation, the concentration of nitrate+nitrite (NOₓ) increased in all layers of the sediment. The concentrations of ammonium, phosphate and silicate, in contrast, decreased during the incubation
<table><caption>**Table 1.** Core-wise nutrient fluxes used as response variables, and macrofauna densities and biomasses, and bioturbation parameters used as predictors in the DistLM-analysis. Note that the values for *Marenzelleria* include both the natural community and the worms added as part of the density treatments. Bioturbation parameters $D_b^N$ = biodiffusion coefficient, SR = % of surface reworked, MPD = maximum penetration depth, BI = bioirrigation.</caption><thead><tr><th rowspan="2">Core</th><th>NO<sub>x</sub></th><th>NH<sub>4</sub><sup>+</sup></th><th>PO<sub>4</sub><sup>3-</sup></th><th>Si<sup>4+</sup></th><th>M. arctica</th><th>M. neg+vir</th><th>Marenzelleria</th><th>C. glaucum</th><th>M. balthica</th><th>H. diversicolo</th><th>Hydrobiidae</th><th>D<sub>b</sub><sup>N</sup></th><th>SR</th><th>MPD</th><th>BI</th></tr><tr><th colspan="2">mmol m<sup>-2</sup> d<sup>-1</sup></th><th colspan="2">ind m<sup>-2</sup></th><th colspan="2">ww g m<sup>-2</sup></th><th colspan="2">cm<sup>2</sup> yr<sup>-1</sup></th><th colspan="2">cm<sup>2</sup> yr<sup>-1</sup></th><th colspan="2">cm<sup>2</sup> yr<sup>-1</sup></th><th colspan="2">ml d<sup>-1</sup></th></tr></thead><tbody><tr><td>Ca</td><td>-0.02</td><td>-0.15</td><td>-2.08</td><td>-0.08</td><td>0</td><td>1805</td><td>23.43</td><td>0.00</td><td>155.56</td><td>0.00</td><td>3.65</td><td>0.73</td><td>83.23</td><td>6.50</td><td>14.10</td></tr><tr><td>Cb</td><td>0.05</td><td>0.07</td><td>0.50</td><td>0.62</td><td>211</td><td>1053</td><td>10.27</td><td>3.27</td><td>0.14</td><td>0.00</td><td>7.55</td><td>0.00</td><td>76.83</td><td>5.50</td><td>32.74</td></tr><tr><td>Cc</td><td>0.10</td><td>0.38</td><td>3.37</td><td>0.85</td><td>0</td><td>1625</td><td>15.00</td><td>7.38</td><td>24.62</td><td>0.09</td><td>5.13</td><td>5.84</td><td>82.79</td><td>5.50</td><td>41.35</td></tr><tr><td>Cd</td><td>0.45</td><td>0.46</td><td>86.82</td><td>4.90</td><td>181</td><td>0</td><td>1.93</td><td>0.00</td><td>186.82</td><td>49.69</td><td>5.32</td><td>2.19</td><td>57.56</td><td>8.50</td><td>10.98</td></tr><tr><td>1a</td><td>-0.01</td><td>0.46</td><td>52.85</td><td>0.29</td><td>812</td><td>812</td><td>21.34</td><td>0.00</td><td>238.50</td><td>3.61</td><td>3.19</td><td>5.84</td><td>68.57</td><td>10.00</td><td>78.89</td></tr><tr><td>1b</td><td>0.03</td><td>0.25</td><td>18.64</td><td>2.98</td><td>0</td><td>1264</td><td>1.75</td><td>0.00</td><td>6.84</td><td>4.89</td><td>5.05</td><td>1.46</td><td>74.49</td><td>5.50</td><td>36.31</td></tr><tr><td>1c</td><td>0.07</td><td>0.39</td><td>-4.22</td><td>0.10</td><td>542</td><td>0</td><td>2.71</td><td>9.35</td><td>0.00</td><td>3.77</td><td>5.07</td><td>1.10</td><td>58.86</td><td>8.50</td><td>21.72</td></tr><tr><td>1d</td><td>-0.02</td><td>0.81</td><td>4.59</td><td>-0.86</td><td>650</td><td>975</td><td>28.03</td><td>27.20</td><td>0.00</td><td>28.23</td><td>5.34</td><td>3.29</td><td>77.01</td><td>6.50</td><td>97.63</td></tr><tr><td>2a</td><td>0.09</td><td>0.61</td><td>18.89</td><td>0.78</td><td>421</td><td>843</td><td>41.75</td><td>6.61</td><td>0.00</td><td>0.45</td><td>4.68</td><td>0.73</td><td>79.12</td><td>5.50</td><td>95.07</td></tr><tr><td>2b</td><td>-0.06</td><td>0.14</td><td>8.52</td><td>0.41</td><td>316</td><td>948</td><td>7.17</td><td>2.40</td><td>183.93</td><td>0.00</td><td>2.33</td><td>0.37</td><td>69.50</td><td>5.50</td><td>37.63</td></tr><tr><td>2c</td><td>0.05</td><td>0.25</td><td>-2.31</td><td>0.04</td><td>1011</td><td>1517</td><td>13.32</td><td>2.33</td><td>0.00</td><td>0.00</td><td>2.69</td><td>1.10</td><td>68.69</td><td>4.50</td><td>22.53</td></tr><tr><td>2d</td><td>0.01</td><td>0.85</td><td>14.60</td><td>-0.45</td><td>963</td><td>481</td><td>10.79</td><td>2.24</td><td>0.02</td><td>0.00</td><td>3.09</td><td>1.46</td><td>71.63</td><td>12.00</td><td>57.80</td></tr><tr><td>3a</td><td>0.00</td><td>0.44</td><td>-1.69</td><td>0.17</td><td>1341</td><td>1006</td><td>33.57</td><td>2.17</td><td>0.00</td><td>2.94</td><td>2.94</td><td>2.19</td><td>80.05</td><td>7.50</td><td>14.27</td></tr><tr><td>3b</td><td>0.02</td><td>0.27</td><td>-1.60</td><td>1.26</td><td>963</td><td>481</td><td>4.95</td><td>1.23</td><td>0.02</td><td>0.00</td><td>1.23</td><td>0.37</td><td>65.73</td><td>3.50</td><td>54.24</td></tr><tr><td>3c</td><td>-0.02</td><td>0.51</td><td>0.57</td><td>-0.26</td><td>722</td><td>1805</td><td>35.18</td><td>11.99</td><td>124.48</td><td>0.00</td><td>7.91</td><td>1.46</td><td>69.93</td><td>8.50</td><td>45.77</td></tr><tr><td>4a</td><td>0.08</td><td>0.97</td><td>16.56</td><td>0.32</td><td>2949</td><td>843</td><td>34.01</td><td>8.10</td><td>173.39</td><td>0.00</td><td>5.23</td><td>4.02</td><td>80.06</td><td>5.50</td><td>67.37</td></tr><tr><td>4b</td><td>0.04</td><td>0.89</td><td>1.74</td><td>1.54</td><td>4198</td><td>1399</td><td>25.72</td><td>3.05</td><td>209.76</td><td>0.00</td><td>2.33</td><td>2.19</td><td>86.64</td><td>7.50</td><td>90.18</td></tr><tr><td>4c</td><td>0.02</td><td>0.93</td><td>8.11</td><td>-0.07</td><td>2785</td><td>464</td><td>20.63</td><td>3.92</td><td>0.00</td><td>14.10</td><td>0.99</td><td>1.10</td><td>80.98</td><td>10.00</td><td>60.01</td></tr><tr><td>4d</td><td>-0.02</td><td>1.04</td><td>3.76</td><td>0.93</td><td>2437</td><td>4062</td><td>65.47</td><td>2.92</td><td>0.00</td><td>0.00</td><td>2.40</td><td>4.75</td><td>88.35</td><td>8.50</td><td>131.25</td></tr><tr><td>5b</td><td>-0.09</td><td>1.39</td><td>12.11</td><td>1.65</td><td>9749</td><td>0</td><td>43.86</td><td>4.39</td><td>14.88</td><td>0.00</td><td>4.04</td><td>1.46</td><td>91.34</td><td>4.50</td><td>140.17</td></tr><tr><td>5c</td><td>-0.01</td><td>-0.70</td><td>-8.27</td><td>-2.19</td><td>433</td><td>650</td><td>14.58</td><td>42.64</td><td>9.64</td><td>164.87</td><td>3.63</td><td>2.56</td><td>90.30</td><td>8.50</td><td>103.81</td></tr><tr><td>average</td><td>0.04</td><td>0.49</td><td>11.02</td><td>0.62</td><td>1461</td><td>1049</td><td>21.69</td><td>6.72</td><td>63.27</td><td>12.98</td><td>3.99</td><td>2.11</td><td>74.70</td><td>7.05</td><td>59.71</td></tr><tr><td>min</td><td>-0.09</td><td>-0.70</td><td>-8.27</td><td>-2.19</td><td>0</td><td>0</td><td>1.75</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.99</td><td>0.00</td><td>25.83</td><td>3.50</td><td>10.98</td></tr><tr><td>max</td><td>0.45</td><td>1.39</td><td>86.82</td><td>4.90</td><td>9749</td><td>4062</td><td>65.47</td><td>42.64</td><td>238.50</td><td>164.87</td><td>7.91</td><td>5.84</td><td>91.34</td><td>12.00</td><td>140.17</td></tr></tbody></table>

surviving was between 67 and 100%. The realized densities and biomasses of the most common infauna are presented in Table 1.

**Solute fluxes.** NOₓ fluxes in the experimental cores ranged from -0.092 mmol m⁻² d⁻¹ (high-density core 5b) to 0.446 mmol m⁻² d⁻¹ (unmanipulated core Cd), ammonium fluxes ranged from -0.704 mmol m⁻² d⁻¹ (high-density core 5c) to 1.388 (high-density treatment 5b), phosphate fluxes ranged from -0.00827 mmol m⁻² d⁻¹ (high-density core 5c) to 0.0866 mmol m⁻² d⁻¹ measured in one of the controls (Cd), and silicate fluxes ranged from 2.186 mmol m⁻² d⁻¹ (in 5c) to 4.897 mmol m⁻² d⁻¹ (in Cd). Core-wise solute fluxes are presented in Table 1.

**Bioturbation metrics.** At the end of the experiment, burrow networks created by *Marenzelleria* spp. were seen along the core walls. Luminophore tracers were buried down into these burrows (Fig. 1A). Typical burrows created by *H. diversicolor* and siphonal gallery networks generated by *M. balthica* during deposit feeding and associated with tracer burial were also apparent. We could also observe fecal-pellet strings, typical of *Marenzelleria* spp., around the openings of the burrows at the sediment surface (Fig. 1B). These pellets were composed of both sediment and luminophore tracers, indicating that *Marenzelleria* spp. individuals had deposit-fed at the surface during the experiment.

The biodiffusion coefficient $D_b^N$ ranged from 0 (unmanipulated core Cb) to 5.84 cm² d⁻¹ (low-density manipulated core 1a). The percentage of surface reworked (SR) ranged from 57.6% (low-density manipulated core 1c) to 91.3% (high-density manipulated core 5b). The maximum penetration depth (MPD) ranged from 3.5 cm (core 3b) to 12.0 cm (core 2d). Bioirrigation (BI) had a minimum of 10.98 ml d⁻¹ in the unmanipulated core Cd, and a maximum of 140.17 ml d⁻¹ in the high-density treated core 5b. Core-wise bioturbation metrics are presented in Table 1.

**Pore water profiles and ΔNH₄⁺, ΔNOₓ, ΔPO₄³⁻ and ΔSi⁴⁺ during dark incubation.** At the start of the incubation, the pore water profiles of the control core Ca (Fig. 2a,b) showed profiles with increasing concentrations of ammonium, phosphate and silicate, and decreasing concentration of nitrate+nitrite from the surface down. During the incubation, the concentration of nitrate+nitrite (NOₓ) increased in all layers of the sediment. The concentrations of ammonium, phosphate and silicate, in contrast, decreased during the incubation
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Photographs showing luminophore tracer fate consecutive to *Marenzelleria* spp. bioturbation and feeding activities at the end of the experiment. (A) Occurrence of *Marenzelleria* burrows (indicated by white arrows) burying tracers down in the sediment along glass-wall of core 5b. White line indicates the sediment-water interface (SWI). Sediment oxidized layer (ox) had a light brown color whereas reduced layer (red) was black. (B) Close-up photo from above core 2a showing the sediment surface with luminophore tracers (in pink and yellow), a *Marenzelleria* burrow opening (b) and fecal-pellet strings (p).</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Pore water profiles of the different solutes at the start (a) and end (b) of the dark incubation of the sediment cores Ca, 1a, 2a, 3a, 4a and 5a. The blue color marks the water column. C = no worms added, 1 = 3 worms, 2 = 6 worms, 3 = 12 worms, 4 = 24 worms and 5 = 48 worms added.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Photographs showing luminophore tracer fate consecutive to *Marenzelleria* spp. bioturbation and feeding activities at the end of the experiment. (A) Occurrence of *Marenzelleria* burrows (indicated by white arrows) burying tracers down in the sediment along glass-wall of core 5b. White line indicates the sediment-water interface (SWI). Sediment oxidized layer (ox) had a light brown color whereas reduced layer (red) was black. (B) Close-up photo from above core 2a showing the sediment surface with luminophore tracers (in pink and yellow), a *Marenzelleria* burrow opening (b) and fecal-pellet strings (p).</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Pore water profiles of the different solutes at the start (a) and end (b) of the dark incubation of the sediment cores Ca, 1a, 2a, 3a, 4a and 5a. The blue color marks the water column. C = no worms added, 1 = 3 worms, 2 = 6 worms, 3 = 12 worms, 4 = 24 worms and 5 = 48 worms added.</figcaption></figure>
<table><caption>**Table 2.** DistLM results for all macrofauna densities and biomasses and bioturbation parameters as predictors of the solutes, and all macrofauna as predictors of the different bioturbation measures. $D_b^N$ = biodiffusion coefficient, SR = % of surface reworked, MPD = maximum penetration depth, BI = bioirrigation. Macrofauna includes *M. arctia*, *M. neglecta* and *M. viridis* density, *Marenzelleria* spp. biomass, *M. balthica* biomass, *C. glaucum* biomass and *H. diversicolor* biomass.</caption><thead><tr><th>Available predictors</th><th>Response</th><th>Selection procedure</th><th>Selected predictors</th><th>AIC</th><th>Pseudo-F</th><th>P</th><th>Prop.</th><th>Cum-R<sup>2</sup></th><th>Correlation with dbRDA-axis</th></tr></thead><tbody><tr><td rowspan="2">Macrofauna + bioturbation</td><td rowspan="2">NO<sub>x</sub></td><td rowspan="2">Sequential (forward)</td><td>Bioirrigation</td><td>194.60</td><td>3.88</td><td>0.06</td><td>0.17</td><td>0.17</td><td>—</td></tr><tr><td>Hediste biomass</td><td>193.23</td><td>3.13</td><td>0.10</td><td>0.12</td><td>0.29</td><td>+</td></tr><tr><td rowspan="2">Macrofauna + bioturbation</td><td rowspan="2">NH<sub>4</sub><sup>+</sup></td><td rowspan="2">Sequential (forward)</td><td>M. arctia density</td><td>245.07</td><td>20.95</td><td>0.00</td><td>0.52</td><td>0.52</td><td>—</td></tr><tr><td>Cerastoderma biomass</td><td>-11.03</td><td>8.07</td><td>0.01</td><td>0.30</td><td>0.30</td><td>+</td></tr><tr><td rowspan="2">Macrofauna + bioturbation</td><td rowspan="2">PO<sub>4</sub><sup>3-</sup></td><td rowspan="2">Sequential (forward)</td><td>D<sub>b</sub><sup>N</sup></td><td>-11.76</td><td>2.50</td><td>0.13</td><td>0.09</td><td>0.38</td><td>—</td></tr><tr><td>M. neg+vir density</td><td>-13.53</td><td>3.34</td><td>0.09</td><td>0.10</td><td>0.48</td><td>+</td></tr><tr><td rowspan="3">Macrofauna+bioturbation</td><td rowspan="3">Si<sup>4+</sup></td><td rowspan="3">Sequential (forward)</td><td>Cerastoderma biomass</td><td>297.19</td><td>12.55</td><td>0.00</td><td>0.40</td><td>0.40</td><td>—</td></tr><tr><td>M. neg+vir density</td><td>296.44</td><td>2.51</td><td>0.13</td><td>0.07</td><td>0.47</td><td>—</td></tr><tr><td>MPD</td><td>296.22</td><td>1.90</td><td>0.18</td><td>0.05</td><td>0.52</td><td>—</td></tr><tr><td rowspan="3">Macrofauna</td><td rowspan="3">MPD</td><td rowspan="3">Sequential (forward)</td><td>Hediste biomass</td><td>-8.12</td><td>10.12</td><td>0.01</td><td>0.35</td><td>0.35</td><td>+</td></tr><tr><td>Macoma biomass</td><td>-10.67</td><td>4.35</td><td>0.06</td><td>0.13</td><td>0.47</td><td>+</td></tr><tr><td>Marenzelleria biomass</td><td>-12.46</td><td>3.36</td><td>0.08</td><td>0.09</td><td>0.56</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">SR</td><td rowspan="2">Sequential (forward)</td><td>Marenzelleria biomass</td><td>-106.09</td><td>12.14</td><td>0.002</td><td>0.39</td><td>0.39</td><td>+</td></tr><tr><td>Marenzelleria biomass</td><td>-20.68</td><td>2.92</td><td>0.11</td><td>0.13</td><td>0.13</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">D<sub>b</sub><sup>N</sup></td><td rowspan="2">Sequential (forward)</td><td>Hediste biomass</td><td>-21.50</td><td>2.59</td><td>0.12</td><td>0.11</td><td>0.24</td><td>+</td></tr><tr><td>Macoma biomass</td><td>-22.13</td><td>2.27</td><td>0.15</td><td>0.09</td><td>0.33</td><td>+</td></tr><tr><td rowspan="3">Macrofauna</td><td rowspan="3">Bioirrigation</td><td rowspan="3">Sequential (forward)</td><td>Marenzelleria biomass</td><td>144.07</td><td>14.57</td><td>0.00</td><td>0.43</td><td>0.43</td><td>+</td></tr><tr><td>M. arctia density</td><td>141.86</td><td>4.00</td><td>0.06</td><td>0.10</td><td>0.54</td><td>+</td></tr><tr><td>Cerastoderma biomass</td><td>140.97</td><td>2.51</td><td>0.13</td><td>0.06</td><td>0.60</td><td>+</td></tr></tbody></table>

in all layers of the sediment. Similar trends were apparent in all treatments, except for cores 3a and 5a. The concentration of NOx decreased in the bottom water during the incubation in all treatments except for cores Ca and 1a, whereas the concentration of phosphate in the bottom water increased in all treatments except for cores Ca and 3a. Ammonium in the bottom water decreased in treatments 1a to 5a, and increased in core Ca, whereas the opposite was true for silicate.

Apart from the manipulated densities of *Marenzelleria* spp., the background communities present were highly variable across cores (see above) with high densities and large organisms potentially affecting pore water concentrations. The pore water control core (Ca), had a high biomass of *M. balthica* (155.6 g/m²) and a relatively high density of Hydrobiidae (2527 ind/m²) but very little other fauna. Core 3a had the third highest density of *Marenzelleria* spp. (2347 ind/m²) and a density of 1805 ind/m² hydrobids, but otherwise very little other fauna, whereas core 5a had a large *Mya arenaria* weighing 6.0 g (1087.7 g/m²), the highest density of *Marenzelleria* spp. (9567 ind/m²), and the second-highest biomass of *M. balthica* (184.1 g/m²).

**Effects of macrofauna on bioturbation and solute fluxes.** Detailed results including specified percentages, AIC and P-values and direction of the multiple partial correlations with the dbRDA-axes for individual predictors are presented in Table 2. Effects of macrofauna on bioturbation metrics, and of macrofauna and bioturbation metrics on solute fluxes are considered direct effects, effects of macrofauna through bioturbation on solute fluxes are considered indirect effects.

**Direct effects of macrofauna on bioturbation metrics.** Macrofauna predicted in total 56% of MPD, 38% of SR, 33% of $D_b^N$, and 60% of bioirrigation (Fig. 3a). The significant predictors accounting for variation in MPD were the biomass of *H. diversicolor*, the biomass of *M. balthica* and the biomass of *Marenzelleria* spp., in SR *Marenzelleria* spp. biomass, and in bioirrigation *Marenzelleria* spp. biomass and the density of *M. arctia*. $D_b^N$ was not significantly predicted by any of the macrofauna at the $P \le 0.10$ -level. The combination of all macrofauna had a positive effect on all the bioturbation metrics (Fig. 3a).

**Direct effects of macrofauna and bioturbation on solute fluxes.** The combination of all macrofauna and bioturbation metrics predicted 29% of the NOx-, 52% of the NH₄⁺-, 61% of the PO₄³⁻-, and 62% of the Si⁴⁺-fluxes (Figs 3b, 4). The significant predictors selected accounting for the variation in NOx-fluxes were bioirrigation and *H. diversicolor* biomass, for NH₄⁺-fluxes the density of *M. arctia*, for PO₄³⁻-fluxes *C. glaucum* biomass and the density of *M. neg+vir* and in Si⁴⁺-fluxes *C. glaucum* biomass (Table 2). Bioirrigation and *H. diversicolor* biomass had opposite effects on NOx-fluxes with bioirrigation decreasing and *H. diversicolor* increasing the fluxes. *M. arctia* increased the fluxes of ammonium. The biomass of *C. glaucum* and the density of *M. neg+vir* decreased the fluxes of phosphate, and the biomass of *C. glaucum*, density of *M. neg+vir*, and MPD increased the silicate fluxes (Figs 3b, 4).
<table><caption>**Table 2.** DistLM results for all macrofauna densities and biomasses and bioturbation parameters as predictors of the solutes, and all macrofauna as predictors of the different bioturbation measures. $D_b^N = \text{biodiffusion coefficient}$, $SR = \%$ of surface reworked, $MPD = \text{maximum penetration depth}$, $BI = \text{bioirrigation}$. Macrofauna includes *M. arctia*, *M. neglecta* and *M. viridis* density, *Marenzelleria* spp. biomass, *M. balthica* biomass, *C. glaucum* biomass and *H. diversicolor* biomass.</caption><thead><tr><th>Available predictors</th><th>Response</th><th>Selection procedure</th><th>Selected predictors</th><th>AIC</th><th>Pseudo-F</th><th>P</th><th>Prop.</th><th>Cum-R<sup>2</sup></th><th>Correlation with dbRDA-axis</th></tr></thead><tbody><tr><td rowspan="2">Macrofauna + bioturbation</td><td rowspan="2">NO<sub>x</sub></td><td rowspan="2">Sequential (forward)</td><td>Bioirrigation</td><td>194.60</td><td>3.88</td><td>0.06</td><td>0.17</td><td>0.17</td><td>—</td></tr><tr><td>Hediste biomass</td><td>193.23</td><td>3.13</td><td>0.10</td><td>0.12</td><td>0.29</td><td>+</td></tr><tr><td rowspan="2">Macrofauna + bioturbation</td><td rowspan="2">NH<sub>4</sub><sup>+</sup></td><td rowspan="2">Sequential (forward)</td><td>M. arctia density</td><td>245.07</td><td>20.95</td><td>0.00</td><td>0.52</td><td>0.52</td><td>—</td></tr><tr><td>Cerastoderma biomass</td><td>-11.03</td><td>8.07</td><td>0.01</td><td>0.30</td><td>0.30</td><td>+</td></tr><tr><td rowspan="2">Macrofauna + bioturbation</td><td rowspan="2">PO<sub>4</sub><sup>3-</sup></td><td rowspan="2">Sequential (forward)</td><td>D<sub>b</sub><sup>N</sup></td><td>-11.76</td><td>2.50</td><td>0.13</td><td>0.09</td><td>0.38</td><td>—</td></tr><tr><td>M. neg+vir density</td><td>-13.53</td><td>3.34</td><td>0.09</td><td>0.10</td><td>0.48</td><td>+</td></tr><tr><td rowspan="2">Macrofauna+bioturbation</td><td rowspan="2">Si<sup>4+</sup></td><td rowspan="2">Sequential (forward)</td><td>Cerastoderma biomass</td><td>297.19</td><td>12.55</td><td>0.00</td><td>0.40</td><td>0.40</td><td>—</td></tr><tr><td>M. neg+vir density</td><td>296.44</td><td>2.51</td><td>0.13</td><td>0.07</td><td>0.47</td><td>—</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">MPD</td><td rowspan="2">Sequential (forward)</td><td>MPD</td><td>296.22</td><td>1.90</td><td>0.18</td><td>0.05</td><td>0.52</td><td>—</td></tr><tr><td>Hediste biomass</td><td>-8.12</td><td>10.12</td><td>0.01</td><td>0.35</td><td>0.35</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">SR</td><td rowspan="2">Sequential (forward)</td><td>Macoma biomass</td><td>-10.67</td><td>4.35</td><td>0.06</td><td>0.13</td><td>0.47</td><td>+</td></tr><tr><td>Marenzelleria biomass</td><td>-12.46</td><td>3.36</td><td>0.08</td><td>0.09</td><td>0.56</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">D<sub>b</sub><sup>N</sup></td><td rowspan="2">Sequential (forward)</td><td>Marenzelleria biomass</td><td>-106.09</td><td>12.14</td><td>0.002</td><td>0.39</td><td>0.39</td><td>+</td></tr><tr><td>Marenzelleria biomass</td><td>-20.68</td><td>2.92</td><td>0.11</td><td>0.13</td><td>0.13</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">Bioirrigation</td><td rowspan="2">Sequential (forward)</td><td>Hediste biomass</td><td>-21.50</td><td>2.59</td><td>0.12</td><td>0.11</td><td>0.24</td><td>+</td></tr><tr><td>Macoma biomass</td><td>-22.13</td><td>2.27</td><td>0.15</td><td>0.09</td><td>0.33</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">Bioirrigation</td><td rowspan="2">Sequential (forward)</td><td>Marenzelleria biomass</td><td>144.07</td><td>14.57</td><td>0.00</td><td>0.43</td><td>0.43</td><td>+</td></tr><tr><td>M. arctia density</td><td>141.86</td><td>4.00</td><td>0.06</td><td>0.10</td><td>0.54</td><td>+</td></tr><tr><td rowspan="2">Macrofauna</td><td rowspan="2">Bioirrigation</td><td rowspan="2">Sequential (forward)</td><td>Cerastoderma biomass</td><td>140.97</td><td>2.51</td><td>0.13</td><td>0.06</td><td>0.60</td><td>+</td></tr><tr><td>Cerastoderma biomass</td><td>140.97</td><td>2.51</td><td>0.13</td><td>0.06</td><td>0.60</td><td>+</td></tr></tbody></table>

in all layers of the sediment. Similar trends were apparent in all treatments, except for cores 3a and 5a. The concentration of NOx decreased in the bottom water during the incubation in all treatments except for cores Ca and 1a, whereas the concentration of phosphate in the bottom water increased in all treatments except for cores Ca and 3a. Ammonium in the bottom water decreased in treatments 1a to 5a, and increased in core Ca, whereas the opposite was true for silicate.

Apart from the manipulated densities of *Marenzelleria* spp., the background communities present were highly variable across cores (see above) with high densities and large organisms potentially affecting pore water concentrations. The pore water control core (Ca), had a high biomass of *M. balthica* (155.6 g/m²) and a relatively high density of Hydrobiidae (2527 ind/m²) but very little other fauna. Core 3a had the third highest density of *Marenzelleria* spp. (2347 ind/m²) and a density of 1805 ind/m² hydrobids, but otherwise very little other fauna, whereas core 5a had a large *Mya arenaria* weighing 6.0 g (1087.7 g/m²), the highest density of *Marenzelleria* spp. (9567 ind/m²), and the second-highest biomass of *M. balthica* (184.1 g/m²).

**Effects of macrofauna on bioturbation and solute fluxes.** Detailed results including specified percentages, AIC and P-values and direction of the multiple partial correlations with the dbRDA-axes for individual predictors are presented in Table 2. Effects of macrofauna on bioturbation metrics, and of macrofauna and bioturbation metrics on solute fluxes are considered direct effects, effects of macrofauna through bioturbation on solute fluxes are considered indirect effects.

**Direct effects of macrofauna on bioturbation metrics.** Macrofauna predicted in total 56% of MPD, 38% of SR, 33% of $D_b^N$, and 60% of bioirrigation (Fig. 3a). The significant predictors accounting for variation in MPD were the biomass of *H. diversicolor*, the biomass of *M. balthica* and the biomass of *Marenzelleria* spp., in SR *Marenzelleria* spp. biomass, and in bioirrigation *Marenzelleria* spp. biomass and the density of *M. arctia*. $D_b^N$ was not significantly predicted by any of the macrofauna at the $P \le 0.10$ -level. The combination of all macrofauna had a positive effect on all the bioturbation metrics (Fig. 3a).

**Direct effects of macrofauna and bioturbation on solute fluxes.** The combination of all macrofauna and bioturbation metrics predicted 29% of the NOx-, 52% of the NH₄⁺-, 61% of the PO₄³⁻-, and 62% of the Si⁴⁺-fluxes (Figs 3b, 4). The significant predictors selected accounting for the variation in NOx-fluxes were bioirrigation and *H. diversicolor* biomass, for NH₄⁺-fluxes the density of *M. arctia*, for PO₄³⁻-fluxes *C. glaucum* biomass and the density of *M. neg+vir* and in Si⁴⁺-fluxes *C. glaucum* biomass (Table 2). Bioirrigation and *H. diversicolor* biomass had opposite effects on NOx-fluxes with bioirrigation decreasing and *H. diversicolor* increasing the fluxes. *M. arctia* increased the fluxes of ammonium. The biomass of *C. glaucum* and the density of *M. neg+vir* decreased the fluxes of phosphate, and the biomass of *C. glaucum*, density of *M. neg+vir*, and MPD increased the silicate fluxes (Figs 3b, 4).
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Direct effects of macrofauna on the different bioturbation metrics (a), and indirect effects of macrofauna through bioturbation on solute fluxes, and direct effects of macrofauna and bioturbation on solute fluxes (b) according to the DistLM models. The numbers indicate the amount of variation accounted for by the set of predictors, the sign indicates the direction of the correlation with the dbRDA -axis.</figcaption></figure>

Discussion

Animal densities, biomasses and functional characteristics as well as their interactions affect sediment processes such as bioturbation and nutrient cycling⁷,¹⁷, and hence any changes in these parameters will influence ecosystem functioning. Changes in the abiotic and biotic environment can change the species-specific expression of functional traits⁴²,⁴³. Increasing biomass of the invasive polychaete *Marenzelleria* spp. significantly enhanced all bioturbation metrics examined, whereas the native fauna only had an effect on one or two individual bioturbation metrics. *Marenzelleria* spp. clearly dominated the bioirrigation pattern by predicting 43 out of 60% of total variation, which seemed to have, in the end, a negligible effect on all other but NOₓ fluxes. Direct density-dependent effects of the taxon on the solute fluxes were demonstrable for NH₄⁺ and PO₄³⁻ and Si⁴⁺ fluxes. This has implications for ecosystem functioning especially in deeper, hypoxia-affected areas, where *Marenzelleria* spp. is sometimes the only macrofaunal taxon present. Particle reworking by all three species has previously been found to be negligible³⁵. Our results, however, show that in the macrofauna community, the biomass of *Marenzelleria* spp. was the only significant predictor of the surface reworked and their biomass was also the predictor accounting for most of the variation in the biodiffusion coefficient, but these predictors were not important in predicting the variation in the solute fluxes. Macrofauna in general enhanced all bioturbation metrics, but individual species affected different bioturbation metrics differently, which in turn affected different solutes. The contrasting effects were associated with the variable importance of the separate bioturbation metrics on the individual solutes. The three different species of *Marenzelleria* had complex effects on the solute fluxes through the contrasting indirect effects mediated by their bioturbation. Future experiments with the three species separated from the beginning would be needed, but since their identification requires genetic analysis, conducting these kinds of experiments in areas with all three species co-occurring is challenging.

Of the individual solute fluxes in this study, phosphate and silicate responded to the same faunal components in opposite ways: the biomass of the cockle *C. glaucum* and the density of *M. neglecta* and *M. viridis* had a direct effect on the fluxes, which in the case of phosphate was a decrease in the efflux, but in the case of silicate an increase in the efflux. Increasing densities of cockles, the most significant predictor for both silicate and phosphate fluxes in this study, tend to increase sediment resuspension and erosion⁴⁴, which may increase e.g. silicate fluxes⁴⁵. Previous laboratory studies have found a positive effect of the bioirrigation of all the three *Marenzelleria*
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Direct effects of macrofauna on the different bioturbation metrics (a), and indirect effects of macrofauna through bioturbation on solute fluxes, and direct effects of macrofauna and bioturbation on solute fluxes (b) according to the DistLM models. The numbers indicate the amount of variation accounted for by the set of predictors, the sign indicates the direction of the correlation with the dbRDA -axis.</figcaption></figure>

Discussion

Animal densities, biomasses and functional characteristics as well as their interactions affect sediment processes such as bioturbation and nutrient cycling⁷,¹⁷, and hence any changes in these parameters will influence ecosystem functioning. Changes in the abiotic and biotic environment can change the species-specific expression of functional traits⁴²,⁴³. Increasing biomass of the invasive polychaete *Marenzelleria* spp. significantly enhanced all bioturbation metrics examined, whereas the native fauna only had an effect on one or two individual bioturbation metrics. *Marenzelleria* spp. clearly dominated the bioirrigation pattern by predicting 43 out of 60% of total variation, which seemed to have, in the end, a negligible effect on all other but NOₓ fluxes. Direct density-dependent effects of the taxon on the solute fluxes were demonstrable for NH₄⁺ and PO₄³⁻ and Si⁴⁺ fluxes. This has implications for ecosystem functioning especially in deeper, hypoxia-affected areas, where *Marenzelleria* spp. is sometimes the only macrofaunal taxon present. Particle reworking by all three species has previously been found to be negligible³⁵. Our results, however, show that in the macrofauna community, the biomass of *Marenzelleria* spp. was the only significant predictor of the surface reworked and their biomass was also the predictor accounting for most of the variation in the biodiffusion coefficient, but these predictors were not important in predicting the variation in the solute fluxes. Macrofauna in general enhanced all bioturbation metrics, but individual species affected different bioturbation metrics differently, which in turn affected different solutes. The contrasting effects were associated with the variable importance of the separate bioturbation metrics on the individual solutes. The three different species of *Marenzelleria* had complex effects on the solute fluxes through the contrasting indirect effects mediated by their bioturbation. Future experiments with the three species separated from the beginning would be needed, but since their identification requires genetic analysis, conducting these kinds of experiments in areas with all three species co-occurring is challenging.

Of the individual solute fluxes in this study, phosphate and silicate responded to the same faunal components in opposite ways: the biomass of the cockle *C. glaucum* and the density of *M. neglecta* and *M. viridis* had a direct effect on the fluxes, which in the case of phosphate was a decrease in the efflux, but in the case of silicate an increase in the efflux. Increasing densities of cockles, the most significant predictor for both silicate and phosphate fluxes in this study, tend to increase sediment resuspension and erosion⁴⁴, which may increase e.g. silicate fluxes⁴⁵. Previous laboratory studies have found a positive effect of the bioirrigation of all the three *Marenzelleria*
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Results from the dbRDA-analysis with macrofauna and bioturbation metrics as predictors of the solute fluxes (see Table 2 for the predictors selected and the signs of the correlations with the dbRDA-axes). The y-axis represents the Euclidean distances between the sampled cores for the selected model.</figcaption></figure>

spp. species on phosphate effluxes, and an even greater enhancement by *M. viridis* and *M. neglecta* compared to *M. arctia*³⁶, which partly contradicts our findings. The combination of *M. viridis* and *M. neglecta* explained more variation than *M. arctia* in phosphate fluxes, hence supporting the findings of Renz and Forster³⁶, but on the other hand bioirrigation was not selected as a significant predictor of phosphate fluxes. However, an increase in the biodiffusion coefficient enhances phosphate efflux, and the biomass of *Marenzelleria* spp. enhances biodiffusion predicting the highest proportion, 13%, of the variation of all macrofauna, hence there might be an enhancement through density-dependent indirect effects on bioturbation. The positive correlation of *C. glaucum* biomass with the dbRDA -axis indicates that their increase decreases the efflux of phosphate. Highest phosphate fluxes tended also to occur when *C. glaucum* was absent, and *H. diversicolor* and *M. balthica* biomasses were high, highlighting the opposite roles of the large bivalves *C. glaucum* and *M. balthica* in this system. Large bivalves are known to have important effects on ecosystem functioning⁴⁶,⁴⁷. The importance of *C. glaucum* living at the surface sediment could be attributed to its role as a link between the benthic and pelagic environments¹⁷. Affecting only the top 2 cm of the sediment³, *C. glaucum* could enhance the buffering capacity of the top-sediment for phosphate by increasing the oxygenation of iron to Fe³⁺ that PO₄³⁻ can adsorb to⁴⁸,⁴⁹, resulting in a decrease in the phosphate efflux. The decrease in phosphate efflux due to increasing densities of *C. glaucum* and *M. neg+vir* might also have been caused by their stimulation of bacterial growth leading to binding of the nutrients in the microbial biomass¹⁹.

Depending e.g. on their bioturbation mode and position in the sediment, different taxa affect the hydrodynamic conditions in the sediment and at the sediment-water -interface differently⁹,⁴⁶,⁵⁰. The functional roles of individual species might change when studied in the presence of other fauna⁴², which could thus change the direction and size of the effect of individual species: porosity of the sediment and the volume of oxygenized sediment can also increase due to activities of bivalves present in the experimental cores¹⁷,⁵¹, possibly modifying the effect of *Marenzelleria* spp. on solute fluxes in the already bioturbated sediment. Joensuu et al.⁵² found that in muddy sediments the densities of *Marenzelleria* spp. correlated negatively with the sediment erosion threshold,
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Results from the dbRDA-analysis with macrofauna and bioturbation metrics as predictors of the solute fluxes (see Table 2 for the predictors selected and the signs of the correlations with the dbRDA-axes). The y-axis represents the Euclidean distances between the sampled cores for the selected model.</figcaption></figure>

spp. species on phosphate effluxes, and an even greater enhancement by *M. viridis* and *M. neglecta* compared to *M. arctia*³⁶, which partly contradicts our findings. The combination of *M. viridis* and *M. neglecta* explained more variation than *M. arctia* in phosphate fluxes, hence supporting the findings of Renz and Forster³⁶, but on the other hand bioirrigation was not selected as a significant predictor of phosphate fluxes. However, an increase in the biodiffusion coefficient enhances phosphate efflux, and the biomass of *Marenzelleria* spp. enhances biodiffusion predicting the highest proportion, 13%, of the variation of all macrofauna, hence there might be an enhancement through density-dependent indirect effects on bioturbation. The positive correlation of *C. glaucum* biomass with the dbRDA -axis indicates that their increase decreases the efflux of phosphate. Highest phosphate fluxes tended also to occur when *C. glaucum* was absent, and *H. diversicolor* and *M. balthica* biomasses were high, highlighting the opposite roles of the large bivalves *C. glaucum* and *M. balthica* in this system. Large bivalves are known to have important effects on ecosystem functioning⁴⁶,⁴⁷. The importance of *C. glaucum* living at the surface sediment could be attributed to its role as a link between the benthic and pelagic environments¹⁷. Affecting only the top 2 cm of the sediment³, *C. glaucum* could enhance the buffering capacity of the top-sediment for phosphate by increasing the oxygenation of iron to Fe³⁺ that PO₄³⁻ can adsorb to⁴⁸,⁴⁹, resulting in a decrease in the phosphate efflux. The decrease in phosphate efflux due to increasing densities of *C. glaucum* and *M. neg+vir* might also have been caused by their stimulation of bacterial growth leading to binding of the nutrients in the microbial biomass¹⁹.

Depending e.g. on their bioturbation mode and position in the sediment, different taxa affect the hydrodynamic conditions in the sediment and at the sediment-water -interface differently⁹,⁴⁶,⁵⁰. The functional roles of individual species might change when studied in the presence of other fauna⁴², which could thus change the direction and size of the effect of individual species: porosity of the sediment and the volume of oxygenized sediment can also increase due to activities of bivalves present in the experimental cores¹⁷,⁵¹, possibly modifying the effect of *Marenzelleria* spp. on solute fluxes in the already bioturbated sediment. Joensuu et al.⁵² found that in muddy sediments the densities of *Marenzelleria* spp. correlated negatively with the sediment erosion threshold,
and positively with the erosion rate. Grazers, like *Peringia ulvae*, also numerous in our cores, can also cause disturbance and loosening of the sediment surface⁵³. Increased resuspension can also lead to increased feeding by the suspension-feeding bivalves⁵⁰, modifying organic matter availability and activity of microbial communities even further. The large role of *C. glaucum* in predicting the phosphate and silicate fluxes, 30 and 40% of variation accounted for, respectively, was nevertheless somewhat surprising given that in some previous experiments they have not been found to greatly affect solute fluxes or the microbial communities as the water pumped through its siphons is not in direct contact with the sediment biogeochemical environment¹⁷. The increase in surface rework- ing with increasing densities of *Marenzelleria* spp. suggests they could indirectly be able to modify the behavior of other species, such as *C. glaucum* or surface deposit feeders. The introduction of *Marenzelleria* spp. could thus have enhanced the strength of benthic-pelagic coupling in these environments through its effects on the native fauna. This is essential for the efficient functioning of the entire ecosystem⁶, and could have implications espe- cially in the undisturbed areas, where the richness of native species is still high.

The relatively low level of variance explained, 29%, for nitrate/nitrite -fluxes by macrofauna and bioturbation compared to the other fluxes, 48 to 52%, indicates that other factors not used as predictors here, such as the microbial community composition⁵⁴, abundance of meiofauna⁵⁵ and availability of ammonium and oxygen⁵⁶, are important for the cycling of these nutrients. Although not quantified here, meiofauna bioturbation is also known to stimulate bacterial denitrification in the Baltic Sea sediments, when large bivalves are absent⁵⁵. Recent observational evidence, however, suggests that the influence of *Marenzelleria* spp. can be variable seasonally, with a potential for a very high importance in predicting NOₓ fluxes at both shallow and deep sites associated with peaking densities, rising temperatures and organic matter input during the spring bloom¹⁰. The most important predictor for NOₓ fluxes, bioirrigation, allows the overlying, oxygenized water to enter the anaerobic sediment zones, where solute exchange processes such as denitrification and nitrification take place, hence its importance in predicting the fluxes. *H. diversicolor* is a gallery-diffuser burrowing deep into the sediment, and is thus capable of irrigating large volumes of sediment. The presence of infauna in general stimulates mineralization, nitrification and denitrification processes and solute transport⁵⁷. The different burrow ventilation and irrigation behaviors of macrofauna lead to differences in the microbial community affecting e.g. nitrification, which is e.g. enhanced in the presence of *H. diversicolor* but not of *M. viridis*⁹. The density of *M. arctia* was also the only significant predictor for ammonium fluxes predicting 52% of the variation, thus it could be the enhancement of the ammo- nium fluxes by increasing densities of this species that also has an indirect effect on the nitrate/nitrite fluxes. *Marenzelleria* spp. biomass and the density of *M. arctia* were also the most important predictors of bioirrigation, suggesting a significant indirect effect on nitrate/nitrite -fluxes through significantly enhancing bioirrigation. In previous experiments, *M. arctia* has been found to enhance ammonium efflux⁵⁸ similarly to our results.

To understand the complex interactions of invasive species with their environment and ambient communities, experiments incorporating natural variability in a range of factors are imperative. This does, however, come with a cost as high variability can make interpretation difficult. In order to meet the assumptions of our statistical models, we often have to remove the most distinctive outliers, and prefer experimental setups with the least number of uncontrolled factors. By doing this we, however, reduce our ability to extrapolate to natural conditions. Even though the average might give us the standard effect, the endpoints of the continuum often contain very interesting information⁵⁹. Due to the uncontrolled community composition, we were forced to remove four replicates from the overall analyses. This points out the difficulty of dealing with many (uncontrolled) sources of variation in these kinds of experiments. Accordingly, future experiments aiming at investigating the density-dependent effects of *Marenzelleria* spp. within the natural community should involve more replicates.

Colonization of previously hypoxic areas by *Marenzelleria* spp. could first lead to increased oxygenation and phosphorus binding⁶⁰ and facilitate the return of the native fauna thus having a major role in the beginning. Later, when the sediment is reoxygenized, the effect of *Marenzelleria* spp. on phosphorus binding is decreased, and their effect on nutrient cycling will be dependent on environmental factors (e.g. organic matter input, temperature¹⁰), and macrofaunal community composition. The combined biomass of all three *Marenzelleria* species affects all the bioturbation metrics, whereas the other dominant members of the community only affect one or a few bioturbation metrics significantly. As a genus they could to some degree compensate for the loss of other members in the community, given their ability to enhance all bioturbation metrics. The three different species, on the other hand, prefer different habitats and differ in their population dynamics and in occupancy of sites spatially and temporally³³, and differ in their effects on the individual bioturbation metrics, indicating that also the effect of the polychaetes on ecosystem functioning will vary accordingly. The uneven distribution of the different *Marenzelleria* species in the cores could have affected the results, however, the proportions of the species densities in general represent those found in nature³³, with highest adult densities expressed by *M. arctia*, and thus the resulting picture is also realistic. We also believe that because changes in *Marenzelleria* densities did not predict all of the variation in all the measured response variables, the range of variation in *Marenzelleria* density was not too high compared to other species present in the cores. The differential direct effects of the three species on the bioturbation metrics, and on the solute fluxes, differ: *M. arctia* seems to have a large effect on nitrogen cycling, whereas *M. neglecta* and *M. viridis* have an effect on phosphorus and silica cycles. Only a few sites will be occupied by all three species at the same time, thus the genus as a whole will only have a large effect on bioturbation at these few sites. Given the spatial separation of the populations in nature, there will also be a spatial separation in their effects on ecosystem functioning, with contrasting effects on different solutes, further depending on the structure of the surrounding community.

According to our results, the species and functional composition of the community affects the expression of the different bioturbation metrics present in the community, leading to changes in the direction and magnitude of the biogeochemical fluxes. Due to the different solutes being affected by different sets of macrofauna and bioturbation metrics, it is important to consider the whole community and not only taxa in isolation in biodiversity-ecosystem functioning experiments, as also interactions between the species and their functions could be significant for the
and positively with the erosion rate. Grazers, like *Peringia ulvae*, also numerous in our cores, can also cause disturbance and loosening of the sediment surface⁵³. Increased resuspension can also lead to increased feeding by the suspension-feeding bivalves⁵⁰, modifying organic matter availability and activity of microbial communities even further. The large role of *C. glaucum* in predicting the phosphate and silicate fluxes, 30 and 40% of variation accounted for, respectively, was nevertheless somewhat surprising given that in some previous experiments they have not been found to greatly affect solute fluxes or the microbial communities as the water pumped through its siphons is not in direct contact with the sediment biogeochemical environment¹⁷. The increase in surface rework- ing with increasing densities of *Marenzelleria* spp. suggests they could indirectly be able to modify the behavior of other species, such as *C. glaucum* or surface deposit feeders. The introduction of *Marenzelleria* spp. could thus have enhanced the strength of benthic-pelagic coupling in these environments through its effects on the native fauna. This is essential for the efficient functioning of the entire ecosystem⁶, and could have implications especially in the undisturbed areas, where the richness of native species is still high.

The relatively low level of variance explained, 29%, for nitrate/nitrite -fluxes by macrofauna and bioturbation compared to the other fluxes, 48 to 52%, indicates that other factors not used as predictors here, such as the microbial community composition⁵⁴, abundance of meiofauna⁵⁵ and availability of ammonium and oxygen⁵⁶, are important for the cycling of these nutrients. Although not quantified here, meiofauna bioturbation is also known to stimulate bacterial denitrification in the Baltic Sea sediments, when large bivalves are absent⁵⁵. Recent observational evidence, however, suggests that the influence of *Marenzelleria* spp. can be variable seasonally, with a potential for a very high importance in predicting NOₓ fluxes at both shallow and deep sites associated with peaking densities, rising temperatures and organic matter input during the spring bloom¹⁰. The most important predictor for NOₓ fluxes, bioirrigation, allows the overlying, oxygenized water to enter the anaerobic sediment zones, where solute exchange processes such as denitrification and nitrification take place, hence its importance in predicting the fluxes. *H. diversicolor* is a gallery-diffuser burrowing deep into the sediment, and is thus capable of irrigating large volumes of sediment. The presence of infauna in general stimulates mineralization, nitrification and denitrification processes and solute transport⁵⁷. The different burrow ventilation and irrigation behaviors of macrofauna lead to differences in the microbial community affecting e.g. nitrification, which is e.g. enhanced in the presence of *H. diversicolor* but not of *M. viridis*⁹. The density of *M. arctia* was also the only significant predictor for ammonium fluxes predicting 52% of the variation, thus it could be the enhancement of the ammonium fluxes by increasing densities of this species that also has an indirect effect on the nitrate/nitrite fluxes. *Marenzelleria* spp. biomass and the density of *M. arctia* were also the most important predictors of bioirrigation, suggesting a significant indirect effect on nitrate/nitrite -fluxes through significantly enhancing bioirrigation. In previous experiments, *M. arctia* has been found to enhance ammonium efflux⁵⁸ similarly to our results.

To understand the complex interactions of invasive species with their environment and ambient communities, experiments incorporating natural variability in a range of factors are imperative. This does, however, come with a cost as high variability can make interpretation difficult. In order to meet the assumptions of our statistical models, we often have to remove the most distinctive outliers, and prefer experimental setups with the least number of uncontrolled factors. By doing this we, however, reduce our ability to extrapolate to natural conditions. Even though the average might give us the standard effect, the endpoints of the continuum often contain very interesting information⁵⁹. Due to the uncontrolled community composition, we were forced to remove four replicates from the overall analyses. This points out the difficulty of dealing with many (uncontrolled) sources of variation in these kinds of experiments. Accordingly, future experiments aiming at investigating the density-dependent effects of *Marenzelleria* spp. within the natural community should involve more replicates.

Colonization of previously hypoxic areas by *Marenzelleria* spp. could first lead to increased oxygenation and phosphorus binding⁶⁰ and facilitate the return of the native fauna thus having a major role in the beginning. Later, when the sediment is reoxygenized, the effect of *Marenzelleria* spp. on phosphorus binding is decreased, and their effect on nutrient cycling will be dependent on environmental factors (e.g. organic matter input, temperature¹⁰), and macrofaunal community composition. The combined biomass of all three *Marenzelleria* species affects all the bioturbation metrics, whereas the other dominant members of the community only affect one or a few bioturbation metrics significantly. As a genus they could to some degree compensate for the loss of other members in the community, given their ability to enhance all bioturbation metrics. The three different species, on the other hand, prefer different habitats and differ in their population dynamics and in occupancy of sites spatially and temporally³³, and differ in their effects on the individual bioturbation metrics, indicating that also the effect of the polychaetes on ecosystem functioning will vary accordingly. The uneven distribution of the different *Marenzelleria* species in the cores could have affected the results, however, the proportions of the species densities in general represent those found in nature³³, with highest adult densities expressed by *M. arctia*, and thus the resulting picture is also realistic. We also believe that because changes in *Marenzelleria* densities did not predict all of the variation in all the measured response variables, the range of variation in *Marenzelleria* density was not too high compared to other species present in the cores. The differential direct effects of the three species on the bioturbation metrics, and on the solute fluxes, differ: *M. arctia* seems to have a large effect on nitrogen cycling, whereas *M. neglecta* and *M. viridis* have an effect on phosphorus and silica cycles. Only a few sites will be occupied by all three species at the same time, thus the genus as a whole will only have a large effect on bioturbation at these few sites. Given the spatial separation of the populations in nature, there will also be a spatial separation in their effects on ecosystem functioning, with contrasting effects on different solutes, further depending on the structure of the surrounding community.

According to our results, the species and functional composition of the community affects the expression of the different bioturbation metrics present in the community, leading to changes in the direction and magnitude of the biogeochemical fluxes. Due to the different solutes being affected by different sets of macrofauna and bioturbation metrics, it is important to consider the whole community and not only taxa in isolation in biodiversity-ecosystem functioning experiments, as also interactions between the species and their functions could be significant for the
outcome of the sedimentary processes. The magnitude and direction of the effect of *Marenzelleria* spp. on solute fluxes also seems to be dependent on the level of disturbance in the system, as increasing densities of the invader could e.g. enhance phosphate binding under hypoxic but not normoxic conditions. This suggests that in disturbed areas *Marenzelleria* spp. could act as a driver of change, whereas in undisturbed areas, and when disturbed areas possibly are recolonized by the native fauna and recover, its effect on e.g. biogeochemical cycling shifts towards being less important, and more affected by the native community composition. Thus, in the (deeper) seasonally hypoxic areas the consequences of the invasion of *Marenzelleria* spp. for ecosystem functioning can be positive, whereas in areas not affected by hypoxia they may not have added any significant value to the functional capacity of the community²⁵ although seasonally their contribution to nutrient cycling also in these areas can be substantial¹⁰. Including environmental context, observations of natural history characteristics, behavior and interactions thus provides a way forward regarding our understanding of the complex ecosystem effects of invasive species.

## Material and Methods

The experiment was conducted in a temperature-controlled climate room at Tvärminne Zoological Station, University of Helsinki, in October 2014. The sediment cores used as experimental units and the *Marenzelleria* spp. added to the cores were collected from sites in close proximity to the station. All three species of *Marenzelleria* spp. occur at the sites, and the species identities were verified with genetic analyses (Bastrop *et al.* unpublished data).

**Sediment and worm collection.** In total 24 sediment cores (plexiglass cores, 8.4 cm internal diameter) were obtained from a shallow, muddy site at five meters depth by SCUBA. An additional three sediment cores were obtained for sediment grain size, C/N and organic matter analyses. In the laboratory the cores were placed in a tank containing seawater with a flow-through in the tank and separately in each core. The cores were left to settle for one day before adding the worms.

The worms were collected at three different sites close to and at the sediment collection site using a Van Veen grab at 5 to 10 m depth. The sediment was gently sieved and the intact worms collected from the sieve. The worms were then stored in jars with flow-through containing sediment until they were added to the experimental cores the following day.

**Experimental design.** Each of the four replicate experimental blocks consisted of one control core with no added worms, and five density treatments with 3, 6, 12, 24 and 48 added worms per treatment corresponding to densities observed in nature and used in the study of Norkko and Reed *et al.*³⁸. Hereafter the treatments will be referred to with C for control and numbers 1–5 for low to high density treatments, and letters a, b, c and d will denote the replicates. Worms were added into the cores and let to acclimatize and establish their burrows for 5 days. After worm addition the cores were fitted with nets to prevent escaping of the worms before they had burrowed. After addition of particle tracers at the sediment surface (see below), the experiment was let to run for nine days under flow-through incubation in an immersion tank supplied with natural running sea water, with a light-dark cycle corresponding to the ambient light cycle during that period (11 h/13 h light/dark regime). The temperature corresponded to the ambient temperature of the seawater during the experimental period (approximately 14 °C). Temperature, salinity and oxygen concentration of the incoming water and in the cores were followed daily. After nine days, the cores were fitted with lids provided with manual stirrers and incubated in darkness for four hours, and stirred manually with regular intervals to prevent gradients from forming in the cores. At start and end of the incubation samples of oxygen, NOₓ, NH₄⁺, PO₃³⁻ and Si⁴⁺ were obtained from each core. One core per treatment was also fitted with holes on the side for pore water extraction at 0-1, 1-2, 3-4 and 5-6 cm depth in the sediment. Pore water was extracted both at the start and at the end of the incubation and analyzed for NOₓ, NH₄⁺, PO₃³⁻ and Si⁴⁺. After the experiment the cores were sliced at 0.5, 1, 1.5, 2, 3, 4, 5, 6, 7, 8, 9, 11, 13, 15 cm, subsamples were obtained for particle mixing analyses (see below), the visible fauna was collected from the slices, and the rest of the core was sieved through a 0.5 mm sieve. The worms were stored separately in 98% ethanol to allow genetic analyses, and the rest of the macrofauna samples was stored in 70% ethanol and stained with rose Bengal. The samples were sorted under microscope, macrofauna identified at the lowest taxonomic level possible, counted and weighed. *Marenzelleria* spp. specimens were counted and weighed and sent to the University of Rostock for genetic analyses.

**Sediment particle mixing and bioirrigation.** Sediment particle mixing was assessed through incubation of sediment cores using luminophores as sediment particle tracers⁶¹. After the 5-day acclimatization following the worm addition, the flow through was stopped and 2 g DW (Dry Weight) of luminophores (eco-trace®, environmental tracing systems, density = 2.5 g cm⁻³) were suspended, homogenized in seawater and spread at the sediment surface using a Pasteur pipette. Two size fractions of luminophores were used (“mud” with particle diameter between 10 and 70 µm and “sand” between 125 and 250 µm) in a ratio of 1.5 g “mud” and 0.5 g “sand” mimicking the grain size at the sediment collection site. Luminophores were allowed to settle for 1 h before flow-through was restarted. The incubation lasted 9 days.

At the end of incubation, a photograph of the sediment surface from above was taken. From this, the percentage of surface reworked (SR) was obtained by subtracting the surface still occupied by luminophores from the total surface area using image analysis and dividing this by the total surface area (see below). Cores were subsequently sliced (see above), slices were homogenized and an approx. 30 g aliquot of sediment sampled for luminophore counting after ensuring that no macrofauna were trapped. 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 corresponding to a single slice using image analysis software⁶². The relative concentrations of luminophores in each slice were then used to compute corresponding vertical depth profiles. These profiles were used for: (1) the determination of the Maximum Penetration Depth of
outcome of the sedimentary processes. The magnitude and direction of the effect of *Marenzelleria* spp. on solute fluxes also seems to be dependent on the level of disturbance in the system, as increasing densities of the invader could e.g. enhance phosphate binding under hypoxic but not normoxic conditions. This suggests that in disturbed areas *Marenzelleria* spp. could act as a driver of change, whereas in undisturbed areas, and when disturbed areas possibly are recolonized by the native fauna and recover, its effect on e.g. biogeochemical cycling shifts towards being less important, and more affected by the native community composition. Thus, in the (deeper) seasonally hypoxic areas the consequences of the invasion of *Marenzelleria* spp. for ecosystem functioning can be positive, whereas in areas not affected by hypoxia they may not have added any significant value to the functional capacity of the community²⁵ although seasonally their contribution to nutrient cycling also in these areas can be substantial¹⁰. Including environmental context, observations of natural history characteristics, behavior and interactions thus provides a way forward regarding our understanding of the complex ecosystem effects of invasive species.

## Material and Methods

The experiment was conducted in a temperature-controlled climate room at Tvärminne Zoological Station, University of Helsinki, in October 2014. The sediment cores used as experimental units and the *Marenzelleria* spp. added to the cores were collected from sites in close proximity to the station. All three species of *Marenzelleria* spp. occur at the sites, and the species identities were verified with genetic analyses (Bastrop *et al.* unpublished data).

**Sediment and worm collection.** In total 24 sediment cores (plexiglass cores, 8.4 cm internal diameter) were obtained from a shallow, muddy site at five meters depth by SCUBA. An additional three sediment cores were obtained for sediment grain size, C/N and organic matter analyses. In the laboratory the cores were placed in a tank containing seawater with a flow-through in the tank and separately in each core. The cores were left to settle for one day before adding the worms.

The worms were collected at three different sites close to and at the sediment collection site using a Van Veen grab at 5 to 10 m depth. The sediment was gently sieved and the intact worms collected from the sieve. The worms were then stored in jars with flow-through containing sediment until they were added to the experimental cores the following day.

**Experimental design.** Each of the four replicate experimental blocks consisted of one control core with no added worms, and five density treatments with 3, 6, 12, 24 and 48 added worms per treatment corresponding to densities observed in nature and used in the study of Norkko and Reed *et al.*³⁸. Hereafter the treatments will be referred to with C for control and numbers 1–5 for low to high density treatments, and letters a, b, c and d will denote the replicates. Worms were added into the cores and let to acclimatize and establish their burrows for 5 days. After worm addition the cores were fitted with nets to prevent escaping of the worms before they had buried. After addition of particle tracers at the sediment surface (see below), the experiment was let to run for nine days under flow-through incubation in an immersion tank supplied with natural running sea water, with a light-dark cycle corresponding to the ambient light cycle during that period (11 h/13 h light/dark regime). The temperature corresponded to the ambient temperature of the seawater during the experimental period (approximately 14 °C). Temperature, salinity and oxygen concentration of the incoming water and in the cores were followed daily. After nine days, the cores were fitted with lids provided with manual stirrers and incubated in darkness for four hours, and stirred manually with regular intervals to prevent gradients from forming in the cores. At start and end of the incubation samples of oxygen, NOₓ, NH₄⁺, PO₃³⁻ and Si⁴⁺ were obtained from each core. One core per treatment was also fitted with holes on the side for pore water extraction at 0-1, 1-2, 3-4 and 5-6 cm depth in the sediment. Pore water was extracted both at the start and at the end of the incubation and analyzed for NOₓ, NH₄⁺, PO₃³⁻ and Si⁴⁺. After the experiment the cores were sliced at 0.5, 1, 1.5, 2, 3, 4, 5, 6, 7, 8, 9, 11, 13, 15 cm, subsamples were obtained for particle mixing analyses (see below), the visible fauna was collected from the slices, and the rest of the core was sieved through a 0.5 mm sieve. The worms were stored separately in 98% ethanol to allow genetic analyses, and the rest of the macrofauna samples was stored in 70% ethanol and stained with rose Bengal. The samples were sorted under microscope, macrofauna identified at the lowest taxonomic level possible, counted and weighed. *Marenzelleria* spp. specimens were counted and weighed and sent to the University of Rostock for genetic analyses.

**Sediment particle mixing and bioirrigation.** Sediment particle mixing was assessed through incubation of sediment cores using luminophores as sediment particle tracers⁶¹. After the 5-day acclimatization following the worm addition, the flow through was stopped and 2 g DW (Dry Weight) of luminophores (eco-trace®, environmental tracing systems, density = 2.5 g cm⁻³) were suspended, homogenized in seawater and spread at the sediment surface using a Pasteur pipette. Two size fractions of luminophores were used (“mud” with particle diameter between 10 and 70 µm and “sand” between 125 and 250 µm) in a ratio of 1.5 g “mud” and 0.5 g “sand” mimicking the grain size at the sediment collection site. Luminophores were allowed to settle for 1 h before flow-through was restarted. The incubation lasted 9 days.

At the end of incubation, a photograph of the sediment surface from above was taken. From this, the percentage of surface reworked (SR) was obtained by subtracting the surface still occupied by luminophores from the total surface area using image analysis and dividing this by the total surface area (see below). Cores were subsequently sliced (see above), slices were homogenized and an approx. 30 g aliquot of sediment sampled for luminophore counting after ensuring that no macrofauna were trapped. 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 corresponding to a single slice using image analysis software⁶². The relative concentrations of luminophores in each slice were then used to compute corresponding vertical depth profiles. These profiles were used for: (1) the determination of the Maximum Penetration Depth of
the tracers during the course of the experiment (MPD), and (2) the mathematical fitting of a Continuous Time Random Walk (CTRW) model⁶³ used to derive a single normal biodiffusion coefficient ($D_b N$ in cm² yr⁻¹) value reflecting particle mixing intensity by the resident macrofauna⁶³,⁶⁴.

Bioirrigation was quantified in each core through the measurement of the decrease of an inert bromide (Br⁻) solute tracer spread in the overlying water. A known volume of stock NaBr solution (1 M) was introduced after the incubation to an initial Br⁻- concentration of ca. 10 mM in the overlying water. Incubation then lasted 24 h during which overlying water was stirred using gentle air bubbling. After the addition of NaBr solution 5 ml of overlying water samples were taken at 0 (15 min), 6 and 24 hours. Samples were kept at 4 °C until analysis. Concentration of Br⁻- ions in the water samples was analyzed spectrophotometrically⁶⁵ at Tvärminne Zoological Station (Shimadzu UV-2501 PC) and the relation between bromide concentration in the overlying water and incubation time assessed using simple linear regression⁶⁶,⁶⁷. Bioirrigation rates were then given as a pore water exchange rate Q (in ml.d⁻¹) calculated after⁶⁶.

**Analyses.** *Solute fluxes.* Oxygen samples were analysed by Winkler titration. Water column and pore water solutes were analysed at Tvärminne Zoological Station (Thermo Scientific Aquakem 250). Quantification limits for solutes were $NH_4^+$ 0.0001, $NO_x$ 0.00003, $PO_4^{3-}$ 0.00003 and $Si^{4+}$ 0.0007 mmol l⁻¹. Solute fluxes measured in the experimental cores were calculated from the difference in the concentration between start and end samples as mmol m⁻² d⁻¹. While the oxygen concentration inside the experimental cores was allowed to decrease during the dark incubation following respiratory activities of the benthic community, the oxygen concentration never dropped below 56% saturation, and is therefore unlikely to have affected neither faunal behaviour nor chemical processes.

*Sediment analyses.* Organic matter content of the sediment was determined as loss on ignition (LOI%). The samples were first dried at 60 °C for 48 h and thereafter combusted at 500 °C for 3 h. For grain size analysis, sediment was first placed in 6% hydrogen peroxide -solution to remove the organic matter. Thereafter the sediment was sieved through a series of sieves (2, 1, 0.5, 0.25, 0.125 and 0.063 mm), the different fractions were dried at 60 °C for 48 h or until dry, and the dry weights of the fractions were measured. The C/N -ratio of the sediment was analysed at Tvärminne Zoological Station with a Europa Scientific ANCA-MS 20-20 15 N/13C mass spectrometer after removal of carbonates from the sediment with HCl (aq).

*Statistical analyses.* We used a series of multivariate multiple regression analyses (DISTLM, PERMANOVA + for PRIMER⁶⁸) with distance-based redundancy analysis, dbRDA, to distinguish the effect of *Marenzelleria* spp. and other macrofauna on different solute fluxes and bioturbation parameters. The parameter % of surface reworked was arcsine transformed. The series of models tested the direct effects of bioturbation metrics and macrofauna on solute fluxes and the direct effects of macrofauna on the different bioturbation metrics. Forward and backward selection procedures were used, and Akaike Selection Criterion (AIC) was used as the selection/stopping criterion. The macrofauna used as predictors were the biomass of *M. balthica*, *Cerastoderma glaucum*, *Hediste diversicolor* and *Marenzelleria* spp., and the densities of *M. arctica*, and *M. viridis* + *M. neglecta* + hybrids of these combined (*M. neg+vir* from hereon). Due to highly variable densities or highly correlated density and biomass that would have had a disproportionate effect in the multivariate analysis, we decided to use biomass for large clams and *Hediste* and the combination of all the *Marenzelleria* species. It was also not possible to separate the biomass of *Marenzelleria* on the species level. To explore species-specific density-dependent effects of *Marenzelleria* these were included in the analysis as separate predictors. Despite the high abundance of hydrobio snails (*Peringia ulvae* and *Erobia ventrosa*) in the cores, these did not significantly affect any of the bioturbation metrics or solutes, and were thus removed from the analyses. In the first series of models, the analysis was allowed to select predictors from macrofauna and bioturbation metrics best accounting for variation in the individual solute fluxes. In the second series of models, the analysis was allowed to select predictors from macrofauna best accounting for variation in the different bioturbation metrics. The effect of macrofauna alone excluding bioturbation metrics, and of bioturbation metrics alone excluding macrofauna on solute fluxes was also tested but since the results from these analyses corresponded to the analyses with macrofauna and bioturbation metrics combined, these will not be further discussed. In the analysis, the response variable is a matrix based on Euclidean distances between the core-wise samples of either the different solute fluxes or the bioturbation metrics. The predictor variables are automatically standardized in the analysis, howerer transformation was required for some variables to include non-linear responses. The biomasses of *C. glaucum* and *H. diversicolor* were log(x + 1) -transformed, whereas the densities of *Marenzelleria* spp. and the biomass of *M. balthica* were square root -transformed. The biomass of *Marenzelleria* spp. did not require transformation. The response variable phosphate was log(x + 10) transformed. The other solutes or the bioturbation metrics were not transformed. The relative impact of the predictors was assessed by examining the direction and magnitude of the correlation coefficients. Therefore we believe the significance level 0.1 is sufficient to regard the predictor variables as possibly having an effect on the response variables. The pore water profiles were analysed visually.

**Data availability.** The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.

**References**

1. Berg, P., Rysgaard, S., Funch, P. & Sejr, M. K. Effects of bioturbation on solutes and solids in marine sediments. *Aquatic Microbial Ecology* **26**, 81–94 (2001).

2. Kristensen, E. & Kostka, J. Macrofaunal burrows and irrigation in marine sediment: microbiological and biogeochemical interactions. *Interactions between macro- and microorganisms in marine sediments*, 125–157 (2005).

3. Mermillod-Blondin, F. & Rosenberg, R. Ecosystem engineering: the impact of bioturbation on biogeochemical processes in marine and freshwater benthic habitats. *Aquatic Sciences* **68**, 434–442, https://doi.org/10.1007/s00027-006-0858-x (2006).
the tracers during the course of the experiment (MPD), and (2) the mathematical fitting of a Continuous Time Random Walk (CTRW) model⁶³ used to derive a single normal biodiffusion coefficient ($D_b N$ in cm² yr⁻¹) value reflecting particle mixing intensity by the resident macrofauna⁶³,⁶⁴.

Bioirrigation was quantified in each core through the measurement of the decrease of an inert bromide (Br⁻) solute tracer spread in the overlying water. A known volume of stock NaBr solution (1 M) was introduced after the incubation to an initial Br⁻- concentration of ca. 10 mM in the overlying water. Incubation then lasted 24 h during which overlying water was stirred using gentle air bubbling. After the addition of NaBr solution 5 ml of overlying water samples were taken at 0 (15 min), 6 and 24 hours. Samples were kept at 4 °C until analysis. Concentration of Br⁻ ions in the water samples was analyzed spectrophotometrically⁶⁵ at Tvärminne Zoological Station (Shimadzu UV-2501 PC) and the relation between bromide concentration in the overlying water and incubation time assessed using simple linear regression⁶⁶,⁶⁷. Bioirrigation rates were then given as a pore water exchange rate Q (in ml.d⁻¹) calculated after⁶⁶.

**Analyses.** *Solute fluxes.* Oxygen samples were analysed by Winkler titration. Water column and pore water solutes were analysed at Tvärminne Zoological Station (Thermo Scientific Aquakem 250). Quantification limits for solutes were $NH_4^+$ 0.0001, $NO_x$ 0.00003, $PO_4^{3-}$ 0.00003 and $Si^{4+}$ 0.0007 mmol l⁻¹. Solute fluxes measured in the experimental cores were calculated from the difference in the concentration between start and end samples as mmol m⁻² d⁻¹. While the oxygen concentration inside the experimental cores was allowed to decrease during the dark incubation following respiratory activities of the benthic community, the oxygen concentration never dropped below 56% saturation, and is therefore unlikely to have affected neither faunal behaviour nor chemical processes.

*Sediment analyses.* Organic matter content of the sediment was determined as loss on ignition (LOI%). The samples were first dried at 60 °C for 48 h and thereafter combusted at 500 °C for 3 h. For grain size analysis, sediment was first placed in 6% hydrogen peroxide -solution to remove the organic matter. Thereafter the sediment was sieved through a series of sieves (2, 1, 0.5, 0.25, 0.125 and 0.063 mm), the different fractions were dried at 60 °C for 48 h or until dry, and the dry weights of the fractions were measured. The C/N -ratio of the sediment was analysed at Tvärminne Zoological Station with a Europa Scientific ANCA-MS 20-20 15 N/13C mass spectrometer after removal of carbonates from the sediment with HCl (aq).

*Statistical analyses.* We used a series of multivariate multiple regression analyses (DISTLM, PERMANOVA + for PRIMER⁶⁸) with distance-based redundancy analysis, dbRDA, to distinguish the effect of *Marenzelleria* spp. and other macrofauna on different solute fluxes and bioturbation parameters. The parameter % of surface reworked was arcsine transformed. The series of models tested the direct effects of bioturbation metrics and macrofauna on solute fluxes and the direct effects of macrofauna on the different bioturbation metrics. Forward and backward selection procedures were used, and Akaike Selection Criterion (AIC) was used as the selection/stopping criterion. The macrofauna used as predictors were the biomass of *M. balthica*, *Cerastoderma glaucum*, *Hediste diversicolor* and *Marenzelleria* spp., and the densities of *M. arctia*, and *M. viridis* + *M. neglecta* + hybrids of these combined (*M. neg+vir* from hereon). Due to highly variable densities or highly correlated density and biomass that would have had a disproportionate effect in the multivariate analysis, we decided to use biomass for large clams and *Hediste* and the combination of all the *Marenzelleria* species. It was also not possible to separate the biomass of *Marenzelleria* on the species level. To explore species-specific density-dependent effects of *Marenzelleria* these were included in the analysis as separate predictors. Despite the high abundance of hydrobio snails (*Peringia ulvae* and *Erobia ventrosa*) in the cores, these did not significantly affect any of the bioturbation metrics or solutes, and were thus removed from the analyses. In the first series of models, the analysis was allowed to select predictors from macrofauna and bioturbation metrics best accounting for variation in the individual solute fluxes. In the second series of models, the analysis was allowed to select predictors from macrofauna best accounting for variation in the different bioturbation metrics. The effect of macrofauna alone excluding bioturbation metrics, and of bioturbation metrics alone excluding macrofauna on solute fluxes was also tested but since the results from these analyses corresponded to the analyses with macrofauna and bioturbation metrics combined, these will not be further discussed. In the analysis, the response variable is a matrix based on Euclidean distances between the core-wise samples of either the different solute fluxes or the bioturbation metrics. The predictor variables are automatically standardized in the analysis, howerer transformation was required for some variables to include non-linear responses. The biomasses of *C. glaucum* and *H. diversicolor* were log(x + 1) -transformed, whereas the densities of *Marenzelleria* spp. and the biomass of *M. balthica* were square root -transformed. The biomass of *Marenzelleria* spp. did not require transformation. The response variable phosphate was log(x + 10) transformed. The other solutes or the bioturbation metrics were not transformed. The relative impact of the predictors was assessed by examining the direction and magnitude of the correlation coefficients. Therefore we believe the significance level 0.1 is sufficient to regard the predictor variables as possibly having an effect on the response variables. The pore water profiles were analysed visually.

**Data availability.** The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.

**References**

1. Berg, P., Rysgaard, S., Funch, P. & Sejr, M. K. Effects of bioturbation on solutes and solids in marine sediments. *Aquatic Microbial Ecology* **26**, 81–94 (2001).

2. Kristensen, E. & Kostka, J. Macrofaunal burrows and irrigation in marine sediment: microbiological and biogeochemical interactions. *Interactions between macro-and microorganisms in marine sediments*, 125–157 (2005).

3. Mermillod-Blondin, F. & Rosenberg, R. Ecosystem engineering: the impact of bioturbation on biogeochemical processes in marine and freshwater benthic habitats. *Aquatic Sciences* **68**, 434–442, https://doi.org/10.1007/s00027-006-0858-x (2006).
4. Cloern, J. E. Does the benthos control phytoplankton biomass in South San Francisco Bay? *Marine Ecology Progress Series*, 191–202 (1982).

5. Chauvaud, L., Jean, F., Ragueneau, O. & Thouzeau, G. Long-term variation of the Bay of Brest ecosystem: benthic-pelagic coupling revisited. *Marine Ecology Progress Series*, 35–48 (2000).

6. Griffiths, J. R. et al. The importance of benthic-pelagic coupling for marine ecosystem functioning in a changing world. *Global Change Biology* **23**, 2179–2196 (2017).

7. Braeckman, U. et al. Role of macrofauna functional traits and density in biogeochemical fluxes and bioturbation. *Marine Ecology Progress Series* **399** (2010).

8. Queirós, A. M. et al. Can benthic community structure be used to predict the process of bioturbation in real ecosystems? *Progress in Oceanography* **137**, 559–569 (2015).

9. Vasquez-Cardenas, D., Quintana, C. O., Meysman, F. J., Kristensen, E. & Boschker, H. T. Species-specific effects of two bioturbating polychaetes on sediment chemoautotrophic bacteria. *Marine Ecology Progress Series* **549**, 55–68 (2016).

10. Kauppi, L., Norkko, J., Ikonen, J. & Norkko, A. Seasonal variability in ecosystem functions: quantifying the contribution of invasive species to nutrient cycling in coastal ecosystems. *Marine Ecology Progress Series* **572**, 193–207 (2017).

11. Bertics, V. J. & Ziebis, W. Bioturbation and the role of microniches for sulfate reduction in coastal marine sediments. *Environmental microbiology* **12**, 3022–3034 (2010).

12. Kristensen, E. Organic matter diagenesis at the oxic/anoxic interface in coastal marine sediments, with emphasis on the role of burrowing animals. *Hydrobiologia* **426**, 1–24 (2000).

13. Karlson, K., Bonsdorff, E. & Rosenberg, R. The impact of benthic macrofauna for nutrient fluxes from Baltic Sea sediments. *AMBIO: A Journal of the Human Environment* **36**, 161–167 (2007).

14. Kristensen, E. et al. What is bioturbation? The need for a precise definition for fauna in aquatic sciences. *Marine Ecology Progress Series* **446**, 285–302 (2011).

15. Danovaro, R. et al. Exponential decline of deep-sea ecosystem functioning linked to benthic biodiversity loss. *Current Biology* **18**, 1–8 (2008).

16. Caliman, A., Carneiro, L. S., Bozelli, R. L., Farjalla, V. F. & Esteves, F. A. Bioturbating space enhances the effects of non-additive interactions among benthic ecosystem engineers on cross-habitat nutrient regeneration. *Oikos* **120**, 1639–1648 (2011).

17. Mermillod-Blondin, F., Rosenberg, R., François-Carcaillet, F., Norling, K. & Mauclaire, L. Influence of bioturbation by three benthic infaunal species on microbial communities and biogeochemical processes in marine sediment. *Aquatic Microbial Ecology* **36**, 271–284 (2004).

18. Bremner, J. Species' traits and ecological functioning in marine conservation and management. *Journal of Experimental Marine Biology and Ecology* **366**, 37–47 (2008).

19. Norling, K., Rosenberg, R., Hulth, S., Grémare, A. & Bonsdorff, E. Importance of functional biodiversity and species-specific traits of benthic fauna for ecosystem functions in marine sediment. *Marine Ecology Progress Series* **332**, 11–23 (2007).

20. Belley, R. & Snelgrove, P. V. Relative contributions of biodiversity and environment to benthic ecosystem functioning. *Frontiers in Marine Science* **3**, 242 (2016).

21. Cardinale, B. J., Palmer, M. A. & Collins, S. L. Species diversity enhances ecosystem functioning through interspecific facilitation. *Nature* **415**, 426–429 (2002).

22. Gamfeldt, L. et al. Marine biodiversity and ecosystem functioning: what's known and what's next? *Oikos* **124**, 252–265 (2015).

23. Gurevitch, J. & Padilla, D. K. Are invasive species a major cause of extinctions? *Trends in Ecology & Evolution* **19**, 470–474 (2004).

24. Karlson, A. L., Näslund, J., Rydén, S. & Elmgren, R. Polychaete invader enhances resource utilization in a species-poor system. *Oecologia* **166**, 1055–1065, https://doi.org/10.1007/s00442-011-1936-x (2011).

25. Hewitt, J. E., Norkko, J., Kauppi, L., Villnäs, A. & Norkko, A. Species and functional trait turnover in response to broad-scale change and an invasive species. *Ecosphere* **7** (2016).

26. Eisenhauer, N. et al. Biodiversity–ecosystem function experiments reveal the mechanisms underlying the consequences of biodiversity change in real world ecosystems. *Journal of Vegetation Science* **27**, 1061–1070 (2016).

27. Lohrer, A. M., Thrush, S. F., Hewitt, J. E. & Kraan, C. The up-scaling of ecosystem functions in a heterogeneous world. *Scientific reports* **5** (2015).

28. Norkko, J. et al. Seafloor Ecosystem Function Relationships: *In Situ* Patterns of Change Across Gradients of Increasing Hypoxic Stress. *Ecosystems* **18**, 1424–1439 (2015).

29. Gammal, J., Norkko, J., Pilditch, C. A. & Norkko, A. Coastal Hypoxia and the Importance of Benthic Macrofauna Communities for Ecosystem Functioning. *Estuaries and Coasts* **40**, 457–468 (2017).

30. Kauppi, L., Norkko, A. & Norkko, J. Large-scale species invasion into a low-diversity system: spatial and temporal distribution of the invasive polychaetes Marenzelleria spp. in the Baltic Sea. *Biol Invasions, 1–20* (2015).

31. Bastrop, R. & Blank, M. Multiple Invasions – A Polychaete Genus Enters the Baltic Sea. *Biol Invasions* **8**, 1195–1200, https://doi.org/10.1007/s10530-005-6186-6 (2006).

32. Blank, M., Laine, A. O., Jürss, K. & Bastrop, R. Molecular identification key based on PCR/RFLP for three polychaete sibling species of the genus *Marenzelleria*, and the species' current distribution in the Baltic Sea. *Helgol Mar Res* **62**, 129–141, https://doi.org/10.1007/s10152-007-0081-8 (2008).

33. Kauppi, L., Norkko, A. & Norkko, J. Seasonal population dynamics of the invasive polychaete genus Marenzelleria spp. in contrasting soft-sediment habitats. *Journal of Sea Research* **131**, 46–60 (2018).

34. Dauer, D. M., Maybury, C. A. & Ewing, R. M. Feeding behavior and general ecology of several spionid polychaetes from the Chesapeake Bay. *Journal of Experimental Marine Biology and Ecology* **54**, 21–38 (1981).

35. Renz, J. R. & Forster, S. Are similar worms different? A comparative tracer study on bioturbation in the three sibling species Marenzelleria arctia, M. viridis, and M. neglecta from the Baltic Sea. *Linnol. Oceanogr* **58**, 2046–2058 (2013).

36. Renz, J. R. & Forster, S. Effects of bioirrigation by the three sibling species of Marenzelleria spp. on solute fluxes and porewater nutrient profiles. *Marine Ecology Progress Series* **505**, 145–159 (2014).

37. Jensen, H. S. & Thamdrup, B. In *Proceedings of the Third International Workshop on Phosphorus in Sediments*. 47–59 (Springer).

38. Norkko, J. et al. A welcome can of worms? Hypoxia mitigation by an invasive species. *Global Change Biology* **18**, 422–434, https://doi.org/10.1111/j.1365-2486.2011.02513.x (2012).

39. Vahtera, E. et al. Internal ecosystem feedbacks enhance nitrogen-fixing cyanobacteria blooms and complicate management in the Baltic Sea. *AMBIO: A journal of the Human Environment* **36**, 186–194 (2007).

40. Michaud, E., Desrosiers, G., Mermillod-Blondin, F., Sundby, B. & Stora, G. The functional group approach to bioturbation: the effects of biodiffusers and gallery-diffusers of the Macoma balthica community on sediment oxygen uptake. *Journal of Experimental Marine Biology and Ecology* **326**, 77–88 (2005).

41. Caliman, A. et al. Functional bioturbator diversity enhances benthic-pelagic processes and properties in experimental microcosms. *Journal of the North American Benthological Society* **26**, 450–459 (2007).

42. Wohlgemuth, D., Solan, M. & Godbold, J. A. In *Proc. R. Soc. B. 20162805* (The Royal Society).

43. Wohlgemuth, D., Solan, M. & Godbold, J. A. Specific arrangements of species dominance can be more influential than evenness in maintaining ecosystem process and function. *Scientific reports* **6**, 39325 (2016).

44. Ciutat, A., Widdows, J. & Pope, N. D. Effect of Cerastoderma edule density on near-bed hydrodynamics and stability of cohesive muddy sediments. *Journal of Experimental Marine Biology and Ecology* **346**, 114–126 (2007).
4. Cloern, J. E. Does the benthos control phytoplankton biomass in South San Francisco Bay? *Marine Ecology Progress Series*, 191–202 (1982).

5. Chauvaud, L., Jean, F., Ragueneau, O. & Thouzeau, G. Long-term variation of the Bay of Brest ecosystem: benthic-pelagic coupling revisited. *Marine Ecology Progress Series*, 35–48 (2000).

6. Griffiths, J. R. et al. The importance of benthic-pelagic coupling for marine ecosystem functioning in a changing world. *Global Change Biology* **23**, 2179–2196 (2017).

7. Braeckman, U. et al. Role of macrofauna functional traits and density in biogeochemical fluxes and bioturbation. *Marine Ecology Progress Series* **399** (2010).

8. Queirós, A. M. et al. Can benthic community structure be used to predict the process of bioturbation in real ecosystems? *Progress in Oceanography* **137**, 559–569 (2015).

9. Vasquez-Cardenas, D., Quintana, C. O., Meysman, F. J., Kristensen, E. & Boschker, H. T. Species-specific effects of two bioturbating polychaetes on sediment chemoautotrophic bacteria. *Marine Ecology Progress Series* **549**, 55–68 (2016).

10. Kauppi, L., Norkko, J., Ikonen, J. & Norkko, A. Seasonal variability in ecosystem functions: quantifying the contribution of invasive species to nutrient cycling in coastal ecosystems. *Marine Ecology Progress Series* **572**, 193–207 (2017).

11. Bertics, V. J. & Ziebis, W. Bioturbation and the role of microniches for sulfate reduction in coastal marine sediments. *Environmental microbiology* **12**, 3022–3034 (2010).

12. Kristensen, E. Organic matter diagenesis at the oxic/anoxic interface in coastal marine sediments, with emphasis on the role of burrowing animals. *Hydrobiologia* **426**, 1–24 (2000).

13. Karlson, K., Bonsdorff, E. & Rosenberg, R. The impact of benthic macrofauna for nutrient fluxes from Baltic Sea sediments. *AMBIO: A Journal of the Human Environment* **36**, 161–167 (2007).

14. Kristensen, E. et al. What is bioturbation? The need for a precise definition for fauna in aquatic sciences. *Marine Ecology Progress Series* **446**, 285–302 (2011).

15. Danovaro, R. et al. Exponential decline of deep-sea ecosystem functioning linked to benthic biodiversity loss. *Current Biology* **18**, 1–8 (2008).

16. Caliman, A., Carneiro, L. S., Bozelli, R. L., Farjalla, V. F. & Esteves, F. A. Bioturbating space enhances the effects of non-additive interactions among benthic ecosystem engineers on cross-habitat nutrient regeneration. *Oikos* **120**, 1639–1648 (2011).

17. Mermillod-Blondin, F., Rosenberg, R., François-Carcaillet, F., Norling, K. & Mauclaire, L. Influence of bioturbation by three benthic infaunal species on microbial communities and biogeochemical processes in marine sediment. *Aquatic Microbial Ecology* **36**, 271–284 (2004).

18. Bremner, J. Species' traits and ecological functioning in marine conservation and management. *Journal of Experimental Marine Biology and Ecology* **366**, 37–47 (2008).

19. Norling, K., Rosenberg, R., Hulth, S., Grémare, A. & Bonsdorff, E. Importance of functional biodiversity and species-specific traits of benthic fauna for ecosystem functions in marine sediment. *Marine Ecology Progress Series* **332**, 11–23 (2007).

20. Belley, R. & Snelgrove, P. V. Relative contributions of biodiversity and environment to benthic ecosystem functioning. *Frontiers in Marine Science* **3**, 242 (2016).

21. Cardinale, B. J., Palmer, M. A. & Collins, S. L. Species diversity enhances ecosystem functioning through interspecific facilitation. *Nature* **415**, 426–429 (2002).

22. Gamfeldt, L. et al. Marine biodiversity and ecosystem functioning: what's known and what's next? *Oikos* **124**, 252–265 (2015).

23. Gurevitch, J. & Padilla, D. K. Are invasive species a major cause of extinctions? *Trends in Ecology & Evolution* **19**, 470–474 (2004).

24. Karlson, A. L., Näslund, J., Rydén, S. & Elmgren, R. Polychaete invader enhances resource utilization in a species-poor system. *Oecologia* **166**, 1055–1065, https://doi.org/10.1007/s00442-011-1936-x (2011).

25. Hewitt, J. E., Norkko, J., Kauppi, L., Villnäs, A. & Norkko, A. Species and functional trait turnover in response to broad-scale change and an invasive species. *Ecosphere* **7** (2016).

26. Eisenhauer, N. et al. Biodiversity–ecosystem function experiments reveal the mechanisms underlying the consequences of biodiversity change in real world ecosystems. *Journal of Vegetation Science* **27**, 1061–1070 (2016).

27. Lohrer, A. M., Thrush, S. F., Hewitt, J. E. & Kraan, C. The up-scaling of ecosystem functions in a heterogeneous world. *Scientific reports* **5** (2015).

28. Norkko, J. et al. Seafloor Ecosystem Function Relationships: *In Situ* Patterns of Change Across Gradients of Increasing Hypoxic Stress. *Ecosystems* **18**, 1424–1439 (2015).

29. Gammal, J., Norkko, J., Pilditch, C. A. & Norkko, A. Coastal Hypoxia and the Importance of Benthic Macrofauna Communities for Ecosystem Functioning. *Estuaries and Coasts* **40**, 457–468 (2017).

30. Kauppi, L., Norkko, A. & Norkko, J. Large-scale species invasion into a low-diversity system: spatial and temporal distribution of the invasive polychaetes Marenzelleria spp. in the Baltic Sea. *Biol Invasions, 1–20* (2015).

31. Bastrop, R. & Blank, M. Multiple Invasions – A Polychaete Genus Enters the Baltic Sea. *Biol Invasions* **8**, 1195–1200, https://doi.org/10.1007/s10530-005-6186-6 (2006).

32. Blank, M., Laine, A. O., Jürss, K. & Bastrop, R. Molecular identification key based on PCR/RFLP for three polychaete sibling species of the genus *Marenzelleria*, and the species' current distribution in the Baltic Sea. *Helgol Mar Res* **62**, 129–141, https://doi.org/10.1007/s10152-007-0081-8 (2008).

33. Kauppi, L., Norkko, A. & Norkko, J. Seasonal population dynamics of the invasive polychaete genus Marenzelleria spp. in contrasting soft-sediment habitats. *Journal of Sea Research* **131**, 46–60 (2018).

34. Dauer, D. M., Maybury, C. A. & Ewing, R. M. Feeding behavior and general ecology of several spionid polychaetes from the Chesapeake Bay. *Journal of Experimental Marine Biology and Ecology* **54**, 21–38 (1981).

35. Renz, J. R. & Forster, S. Are similar worms different? A comparative tracer study on bioturbation in the three sibling species Marenzelleria arctica, M. viridis, and M. neglecta from the Baltic Sea. *Linnol. Oceanogr* **58**, 2046–2058 (2013).

36. Renz, J. R. & Forster, S. Effects of bioirrigation by the three sibling species of Marenzelleria spp. on solute fluxes and porewater nutrient profiles. *Marine Ecology Progress Series* **505**, 145–159 (2014).

37. Jensen, H. S. & Thamdrup, B. In *Proceedings of the Third International Workshop on Phosphorus in Sediments*. 47–59 (Springer).

38. Norkko, J. et al. A welcome can of worms? Hypoxia mitigation by an invasive species. *Global Change Biology* **18**, 422–434, https://doi.org/10.1111/j.1365-2486.2011.02513.x (2012).

39. Vahtera, E. et al. Internal ecosystem feedbacks enhance nitrogen-fixing cyanobacteria blooms and complicate management in the Baltic Sea. *AMBIO: A journal of the Human Environment* **36**, 186–194 (2007).

40. Michaud, E., Desrosiers, G., Mermillod-Blondin, F., Sundby, B. & Stora, G. The functional group approach to bioturbation: the effects of biodiffusers and gallery-diffusers of the Macoma balthica community on sediment oxygen uptake. *Journal of Experimental Marine Biology and Ecology* **326**, 77–88 (2005).

41. Caliman, A. et al. Functional bioturbator diversity enhances benthic-pelagic processes and properties in experimental microcosms. *Journal of the North American Benthological Society* **26**, 450–459 (2007).

42. Wohlgemuth, D., Solan, M. & Godbold, J. A. In *Proc. R. Soc. B. 20162805* (The Royal Society).

43. Wohlgemuth, D., Solan, M. & Godbold, J. A. Specific arrangements of species dominance can be more influential than evenness in maintaining ecosystem process and function. *Scientific reports* **6**, 39325 (2016).

44. Ciutat, A., Widdows, J. & Pope, N. D. Effect of Cerastoderma edule density on near-bed hydrodynamics and stability of cohesive muddy sediments. *Journal of Experimental Marine Biology and Ecology* **346**, 114–126 (2007).
45. Katz, T. et al. The silica cycle in a Northeast Pacific fjord; the role of biological resuspension. *Progress in Oceanography* **147**, 10–21 (2016).

46. Thrush, S. F., Hewitt, J. E., Gibbs, M., Lundquist, C. & Norkko, A. Functional role of large organisms in intertidal communities: community effects and ecosystem function. *Ecosystems* **9**, 1029–1040 (2006).

47. Norkko, A., Villnäs, A., Norkko, J., Valanko, S. & Pilditch, C. Size matters: implications of the loss of large individuals for ecosystem function. *Scientific reports* **3**, 2646 (2013).

48. Sundbäck, K. & Granéli, W. Influence of microphytobenthos on the nutrient flux between sediment and water: a laboratory study. *Marine Ecology Progress Series*, **63**, 63–69 (1988).

49. Grace, M. R., Scicluna, T. R., Vithana, C. L., Symes, P. & Lansdown, K. P. Biogeochemistry and cyanobacterial blooms: investigating the relationship in a shallow, polymictic, temperate lake. *Environmental Chemistry* **7**, 443–456, https://doi.org/10.1071/EN10042 (2010).

50. Jones, H. F., Pilditch, C. A., Bryan, K. R. & Hamilton, D. P. Effects of infaunal bivalve density and flow speed on clearance rates and near-bed hydrodynamics. *Journal of Experimental Marine Biology and Ecology* **401**, 20–28 (2011).

51. Michaud, E. et al. Spatial interactions in the Macoma balthica community control biogeochemical fluxes at the sediment-water interface and microbial abundances. *Journal of marine research* **67**, 43–70 (2009).

52. Joensuu, M. et al. Sediment properties, biota, and local habitat structure explain variation in the erodibility of coastal sediments. *Limnology and Oceanography* **63**, 173–186 (2018).

53. Orvain, F., Guizien, K., Lefebvre, S., Bréret, M. & Dupuy, C. Relevance of macrozoobenthic grazers to understand the dynamic behaviour of sediment erodibility and microphytobenthos resuspension in sunny summer conditions. *Journal of Sea Research* **92**, 46–55 (2014).

54. Foshtomi, M. Y. et al. The Link between Microbial Diversity and Nitrogen Cycling in Marine Sediments Is Modulated by Macrofaunal Bioturbation. *PloS one* **10**, e0130116 (2015).

55. Bonaglia, S., Nascimento, F. A., Bartoli, M., Klawonn, I. & Brüchert, V. Meiofauna increases bacterial denitrification in marine sediments. *Nature communications* **5**, 5133 (2014).

56. Rysgaard, S., Risgaard-Petersen, N., Niels Peter, S., Kim, J. & Lars Peter, N. Oxygen regulation of nitrification and denitrification in sediments. *Limnology and Oceanography* **39**, 1643–1652 (1994).

57. Kristensen, E. & Blackburn, T. The fate of organic carbon and nitrogen in experimental marine sediment systems: influence of bioturbation and anoxia. *Journal of Marine Research* **45**, 231–257 (1987).

58. Hietanen, S., Laine, A. O. & Lukkari, K. The complex effects of the invasive polychaetes *Marenzelleria* spp. on benthic nutrient dynamics. *Journal of Experimental Marine Biology and Ecology* **352**, 89–102, https://doi.org/10.1016/j.jembe.2007.07.018 (2007).

59. Hewitt, J., Thrush, S. & Dayton, P. Habitat variation, species diversity and ecological functioning in a marine system. *Journal of Experimental Marine Biology and Ecology* **366**, 116–122 (2008).

60. Bonaglia, S. et al. Effect of reoxygenation and *Marenzelleria* spp. bioturbation on Baltic Sea sediment metabolism. *Marine Ecology Progress Series* **482**, 43–55 (2013).

61. Mahaut, M. & Graf, G. A luminophore tracer technique for bioturbation studies. *Oceanologica acta* **10**, 323–328 (1987).

62. Maire, O. et al. Effects of food availability on sediment reworking in Abra ovata and A. nitida. *Marine Ecology Progress Series* **319**, 136–153 (2006).

63. Meysman, F. J., Malyuga, V. S., Boudreau, B. P. & Middelburg, J. J. A generalized stochastic approach to particle dispersal in soils and sediments. *Geochimica et Cosmochimica Acta* **72**, 3460–3478 (2008).

64. Bernard, G. et al. Comparative study of sediment particle mixing in a Zostera noltei meadow and a bare sediment mudflat. *Marine Ecology Progress Series* **514**, 71–86 (2014).

65. Presley, B. & Claypool, G. Techniques for analyzing interstitial water samples. *Part 1*, 1749–1755 (1971).

66. Rao, A. M., Malkin, S. Y., Montserrat, F. & Meysman, F. J. Alkalinity production in intertidal sands intensified by lugworm bioirrigation. *Estuarine, coastal and shelf science* **148**, 36–47 (2014).

67. De Smet, B., Braeckman, U., Soetaert, K., Vincx, M. & Vanaverbeke, J. Predator effects on the feeding and bioirrigation activity of ecosystem-engineered Lanice conchilega reefs. *Journal of Experimental Marine Biology and Ecology* **475**, 31–37 (2016).

68. Anderson, M., Gorley, R. N. & Clarke, R. K. *Permanova+ for Primer: Guide to software and statistical methods*. (Primer-E Ltd, 2008).

## Acknowledgements

This study was supported financially by Jenny and Antti Wihuri Foundation (grant to LK), Walter and Andrée de Nottbeck Foundation (grant to GB), the BONUS COCOA-project supported by BONUS (Art 185), funded jointly by the EU and the Academy of Finland, and the University of Helsinki (3 yr grant to JN). We thank the personnel at Tvärminne Zoological station for their invaluable help in field sampling and laboratory analyses. We also wish to thank Andrew Lohrer for valuable comments on the manuscript and Alicia Romero-Ramirez for technical support on image analyses.

## Author Contributions

L.K. and G.B. designed the study, performed the field sampling, ran the experiment, performed laboratory analyses, analyzed the data and wrote the paper, J.N. and A.N. designed the study and wrote the paper, R.B. performed the genetic analyses and wrote the paper.

## Additional Information

**Competing Interests:** The authors declare no competing interests.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_6.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2018
45. Katz, T. et al. The silica cycle in a Northeast Pacific fjord; the role of biological resuspension. *Progress in Oceanography* **147**, 10–21 (2016).

46. Thrush, S. F., Hewitt, J. E., Gibbs, M., Lundquist, C. & Norkko, A. Functional role of large organisms in intertidal communities: community effects and ecosystem function. *Ecosystems* **9**, 1029–1040 (2006).

47. Norkko, A., Villnäs, A., Norkko, J., Valanko, S. & Pilditch, C. Size matters: implications of the loss of large individuals for ecosystem function. *Scientific reports* **3**, 2646 (2013).

48. Sundbäck, K. & Granéli, W. Influence of microphytobenthos on the nutrient flux between sediment and water: a laboratory study. *Marine Ecology Progress Series*, **63**, 63–69 (1988).

49. Grace, M. R., Scicluna, T. R., Vithana, C. L., Symes, P. & Lansdown, K. P. Biogeochemistry and cyanobacterial blooms: investigating the relationship in a shallow, polymictic, temperate lake. *Environmental Chemistry* **7**, 443–456, https://doi.org/10.1071/EN10042 (2010).

50. Jones, H. F., Pilditch, C. A., Bryan, K. R. & Hamilton, D. P. Effects of infaunal bivalve density and flow speed on clearance rates and near-bed hydrodynamics. *Journal of Experimental Marine Biology and Ecology* **401**, 20–28 (2011).

51. Michaud, E. et al. Spatial interactions in the Macoma balthica community control biogeochemical fluxes at the sediment-water interface and microbial abundances. *Journal of marine research* **67**, 43–70 (2009).

52. Joensuu, M. et al. Sediment properties, biota, and local habitat structure explain variation in the erodibility of coastal sediments. *Limnology and Oceanography* **63**, 173–186 (2018).

53. Orvain, F., Guizien, K., Lefebvre, S., Bréret, M. & Dupuy, C. Relevance of macrozoobenthic grazers to understand the dynamic behaviour of sediment erodibility and microphytobenthos resuspension in sunny summer conditions. *Journal of Sea Research* **92**, 46–55 (2014).

54. Foshtomi, M. Y. et al. The Link between Microbial Diversity and Nitrogen Cycling in Marine Sediments Is Modulated by Macrofaunal Bioturbation. *PloS one* **10**, e0130116 (2015).

55. Bonaglia, S., Nascimento, F. A., Bartoli, M., Klawonn, I. & Brüchert, V. Meiofauna increases bacterial denitrification in marine sediments. *Nature communications* **5**, 5133 (2014).

56. Rysgaard, S., Risgaard-Petersen, N., Niels Peter, S., Kim, J. & Lars Peter, N. Oxygen regulation of nitrification and denitrification in sediments. *Limnology and Oceanography* **39**, 1643–1652 (1994).

57. Kristensen, E. & Blackburn, T. The fate of organic carbon and nitrogen in experimental marine sediment systems: influence of bioturbation and anoxia. *Journal of Marine Research* **45**, 231–257 (1987).

58. Hietanen, S., Laine, A. O. & Lukkari, K. The complex effects of the invasive polychaetes *Marenzelleria* spp. on benthic nutrient dynamics. *Journal of Experimental Marine Biology and Ecology* **352**, 89–102, https://doi.org/10.1016/j.jembe.2007.07.018 (2007).

59. Hewitt, J., Thrush, S. & Dayton, P. Habitat variation, species diversity and ecological functioning in a marine system. *Journal of Experimental Marine Biology and Ecology* **366**, 116–122 (2008).

60. Bonaglia, S. et al. Effect of reoxygenation and *Marenzelleria* spp. bioturbation on Baltic Sea sediment metabolism. *Marine Ecology Progress Series* **482**, 43–55 (2013).

61. Mahaut, M. & Graf, G. A luminophore tracer technique for bioturbation studies. *Oceanologica acta* **10**, 323–328 (1987).

62. Maire, O. et al. Effects of food availability on sediment reworking in Abra ovata and A. nitida. *Marine Ecology Progress Series* **319**, 136–153 (2006).

63. Meysman, F. J., Malyuga, V. S., Boudreau, B. P. & Middelburg, J. J. A generalized stochastic approach to particle dispersal in soils and sediments. *Geochimica et Cosmochimica Acta* **72**, 3460–3478 (2008).

64. Bernard, G. et al. Comparative study of sediment particle mixing in a Zostera noltei meadow and a bare sediment mudflat. *Marine Ecology Progress Series* **514**, 71–86 (2014).

65. Presley, B. & Claypool, G. Techniques for analyzing interstitial water samples. *Part 1*, 1749–1755 (1971).

66. Rao, A. M., Malkin, S. Y., Montserrat, F. & Meysman, F. J. Alkalinity production in intertidal sands intensified by lugworm bioirrigation. *Estuarine, coastal and shelf science* **148**, 36–47 (2014).

67. De Smet, B., Braeckman, U., Soetaert, K., Vincx, M. & Vanaverbeke, J. Predator effects on the feeding and bioirrigation activity of ecosystem-engineered Lanice conchilega reefs. *Journal of Experimental Marine Biology and Ecology* **475**, 31–37 (2016).

68. Anderson, M., Gorley, R. N. & Clarke, R. K. *Permanova+ for Primer: Guide to software and statistical methods*. (Primer-E Ltd, 2008).

## Acknowledgements

This study was supported financially by Jenny and Antti Wihuri Foundation (grant to LK), Walter and Andrée de Nottbeck Foundation (grant to GB), the BONUS COCOA-project supported by BONUS (Art 185), funded jointly by the EU and the Academy of Finland, and the University of Helsinki (3 yr grant to JN). We thank the personnel at Tvärminne Zoological station for their invaluable help in field sampling and laboratory analyses. We also wish to thank Andrew Lohrer for valuable comments on the manuscript and Alicia Romero-Ramirez for technical support on image analyses.

## Author Contributions

L.K. and G.B. designed the study, performed the field sampling, ran the experiment, performed laboratory analyses, analyzed the data and wrote the paper, J.N. and A.N. designed the study and wrote the paper, R.B. performed the genetic analyses and wrote the paper.

## Additional Information

**Competing Interests:** The authors declare no competing interests.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_7.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.