Introduction

Bivalves are widely distributed across benthic environments and often represent the dominant component of macroinvertebrate communities, both in terms of density and biomass. They filter organic particles from the water column and transfer undigested organic and inorganic particles to the sediment in the form of faeces and pseudofaeces. Next to their filtering and biodepositional functions, bivalves act as bioturbators. Although systems of burrows built by infaunal bivalves are less complex than those built by other macrobenthic bioturbators, their feeding activities—moving about in the sediment and withdrawing and extending their siphons—can modify the physical structure and chemistry of sediments (Dame, 1993; Jones et al., 1994; Gili & Coma, 1998; Gutierrez et al., 2003). Through these modifications, bivalves affect transport rates in the sediment, facilitate oxygen exchange and penetration into sediments (Levinton, 1995; Dame, 1996) and stimulate microbial metabolism, thus affecting other benthic organisms (Olafsson, 2003 and references therein).

Even though it is known that bivalves alter sediment processes, being among the most conspicuous and well-studied invertebrates in the marine environment, their effect on meiofaunal assemblages is relatively poorly studied. Research to date has focused on the effects of single bivalve species locally important in specific geographical regions, sediment types and/or environments (Olafsson & Elmgren, 1991; Olafsson et al., 1993; Warwick et al., 1997; Austen et al., 1998; Olafsson, 2003; Braeckman et al., 2011a). Therefore, these studies cannot be easily extrapolated to other benthic environments and regions. In addition, the available data are often contrasting (for overview see Olafsson, 2003). For instance, Reise (1983) demonstrated that the clam Macoma balthica (Linnaeus, 1758) stimulated turbellarian and nematode abundances, whilst no effect or even a reduced abundance was observed by other authors undertaking the same comparison (Olafsson et al., 1993, 2005, respectively).

Bivalve species vary considerably in terms of mobility, feeding type and activity, method of biodeposits release, and intensity of sediment disturbance. In this study we investigated the effects of three bivalve species, differing in terms of feeding mode (suspension- vs. deposit-feeding), mobility and mode of faeces deposition, on meiobenthos with emphasis on nematodes.

Cerastoderma glaucum (Bruguière, 1789), M. balthica and Mya arenaria Linnaeus, 1758 are the most common burrowing bivalves of the southern Baltic Sea, often dominating macrobenthic communities (Warzocha, 1994, 1995; Piesik et al., 2009). M. balthica is a facultative deposit feeder taking algae, bacteria and detritus from the sediment surface but is able to switch to a suspension-feeding mode depending on food availability (Ólafsson, 1986). It is found buried in the sediment usually at 2–6 cm depth, with its inhalant siphon extended to the sediment surface and the shorter exhalent siphon terminating below the surface. Both M. arenaria and C. glaucum are suspension feeders taking their food out of the water column but their lifestyles in the sediment differ considerably. The soft-shell clam Mya buries deeply (10–25 cm) and leads a sessile lifestyle extending its single, long siphon (inhalant and exhalant siphons are fused) to the sediment surface. It can have a profound effect on sediment biochemistry due to leakage of water from the shell (Hansen et al., 1996). The cockle Cerastoderma has two short, separate siphons and lives actively in the upper few centimetres (0–3 cm), mixing sediment particles and therefore acting as a biodiffuser. Both suspension feeders eject their faeces and pseudofaeces to the sediment surface, whilst Macoma releases them in the subsurface sediment (Ólafsson, 1986; Olafsson et al., 1993). Bivalve biodeposits are rich in carbon and nitrogen (Kautsky & Evans, 1987) and may organically enrich sediments and promote the development of the microbial community. It is hypothesized here that the vertical distribution of meiobenthos and nematodes, and nematode community composition and structure are affected by the presence and species identity of bivalves.

Materials and methods

Characteristics of the study site

Sediment and bivalves for the experiment were collected from a sheltered site located in the inner part of Puck Bay (near Chałupy, Polish Baltic coast), at 70–80 cm water depth. Sediment at the study site was a well-sorted sand dominated by the medium fraction (0.25–0.5 mm; on average 80%), followed by fine (0.1–0.25 mm; 9%) and coarse sand (0.5–1.0 mm; 9%). Organic carbon and total nitrogen (TN) contents averaged 0.06 and 0.02%, respectively, whilst sediment chlorophyll a and phaeophytin concentrations were 1.8 and 0.7 µg g dry sediment−1, respectively (Urban-Malinga et al., 2014).

Sediment and fauna sampling

The experiment was performed in September 2009. Triplicate sediment cores with a surface area of 10 cm2 were collected randomly at the study site to determine the abundance and structure of the meiobenthic community, hereafter referred to as the field community. Sediment cores were sliced immediately on sampling into seven depth layers: 0–1, 1–2, 2–3, 3–4, 4–5, 5–10 and 10–15 cm. These slices were preserved separately in neutral 4% formaldehyde and were further processed in accordance with methods described below for meiofauna in the experiment. Nematode densities and assemblage structure of the field community served as a field control.

Sediment samples for the experiment were taken with a core tube with a 225 cm2 surface area to a depth of 10 cm. The sediment was immediately sieved over 1-mm mesh in a small amount of ambient seawater to exclude macrofauna but to retain all interstitial biota. This may have also excluded larger meiofauna, and although we acknowledge this is not ideal, it was felt to be the most replicable method of eliminating the macrofaunal component. In the laboratory, the sediment was gently homogenized by hand and put into plexiglass cores, 12.3 cm internal diameter and 33 cm long, to a depth of 15 cm. These sediment cores were then placed in a water bath with a total volume of 1100 dm3, connected to a reservoir of 2400 dm3 equipped with a cooling system and an open-loop seawater pumping system. Water for the system was transported directly from the sea over a distance of 50 m, it was filtered through a 2-mm mesh to remove large organic particles and fauna but to retain natural concentrations of phytoplankton and other organic suspensions. Water from the system was distributed to the cores via plastic tubes to facilitate turnover of overlying water (15 cm deep). Outflow water was recycled. Water in the whole system was replaced with fresh seawater once per week. The sediment was allowed to stabilize for 14 days before experiments commenced: preliminary experiments showed that the vertical pattern of meiofauna distribution as recorded at the study site (with the highest meiofauna concentrations at the sediment surface, 0–3 cm) was re-established after this period (Urban-Malinga et al., 2014). The three bivalve species dominating shallow benthic environments of the Gulf of Gdańsk were selected for the study: M. balthica, M. arenaria and C. glaucum (Warzocha, 1994, 1995; Piesik et al., 2009).

Intact specimens were chosen for the experiment: Cerastoderma (9–14 mm, 0.2–0.7 g w. wt.), Macoma (12–15 mm, 0.2–0.4 g w. wt.) and Mya (15–25 mm, 0.3–1.2 g w. wt.). Specimens were added to the microcosms to obtain mono-specific treatments (hereafter referred to as Macoma, Mya and Cerastoderma) of a total bivalve biomass (2.1 g w. wt. per core) similar to that recorded in the field during the study period (a mean of 170 g w. wt. m−2; the three bivalve species selected for this study are the only bivalves at our study site). Biomass (live wet weight) was estimated by weighing each specimen to the nearest 0.1 g. Wet weight excluded water inside the valve (30% of total bivalve weight, pers. obs.) but included shell weight. Densities of Cerastoderma and Mya introduced to the microcosms (on average 7 and 3 ind. per core corresponding to 610 and 280 ind. m−2, respectively) were within the range of their natural densities recorded in the field (i.e. Cerastoderma: 178–889 ind. m−2, Mya: 15–340 ind. m−2), whilst the densities of Macoma (on average 6 ind. per core corresponding to 470 ind. m−2) exceeded these ranges (30–320 ind. m−2). We have used standardized biomass since the densities appropriate for the experiment (e.g. corresponding to average field densities of Cerastoderma or Macoma) would result in unnaturally high densities and biomass of Mya.

Each treatment with bivalves and the control treatment with no macrofauna were performed with three replicate cores.

The majority of specimens buried into the sediment within half an hour of adding to the mesocosms. It was assumed that specimens which had not buried within this time were dead or damaged, and these were replaced. The microcosms were then incubated at 14°C (the ambient temperature at the time of sampling) for 1 month. During this time each microcosm was monitored to control water temperature and overlying water exchange rates and to remove and replace any dead specimens (with an animal of similar size) appearing on the sediment surface. After 1 month, samples for meiofauna and sediment characteristics were taken.

Total organic carbon and total nitrogen contents

Sediment cores for determination of total organic carbon (TOC) and TN contents were collected from each microcosm with a cut-off syringe with an inner diameter of 1.4 cm. Each core was sliced into the following depth layers: 0–1, 1–2, 2–3, 3–4, 4–5 and 5–10 cm. The sediment samples were then frozen at −20°C before further analysis. TOC and TN contents were measured with a Perkin Elmer CHNS/O analyzer after removal of carbonates and drying at 60°C for 24 h.

Pore water nutrient content

For pore water nutrient analysis, one core with an inner diameter of 3.6 cm was sampled from each microcosm and was sliced into the following depth layers: 0–1, 1–2, 2–3, 3–4, 4–5 and 5–10 cm. Samples were frozen at −20°C and then thawed at room temperature prior to analysis. Due to the small amount of pore water obtained from each sediment slice, replicate slices were pooled together. The pore water was extracted with ultra pure water and filtered on a 0.7-µm glass filter. Ammonium and phosphate concentrations were measured immediately after extraction according to standard methods recommended for the Baltic Monitoring Program (Grasshoff et al., 1983).

Meiofauna

One sediment core with an inner diameter of 3.6 cm and a cross-sectional area of 10 cm2 was sampled from each microcosm for determination of the meiofaunal community composition and structure. Sediment cores were sliced immediately into seven depth layers: 0–1, 1–2, 2–3, 3–4, 4–5, 5–10 and 10–15 cm, and each slice was preserved with a neutral 4% formaldehyde solution. Prior to the analysis, samples were first rinsed over a 1-mm mesh to remove macrofauna. Meiofauna was extracted by re-suspending the sediment and decanting the overlying water 10 times over a 38-μm mesh sieve. The fraction retained on the sieve was preserved with a 4% formaldehyde solution and stained with rose bengal. Meiofauna were counted and identified to higher taxon level under a stereomicroscope, and the first 120 nematodes encountered in each sediment slice were picked out and mounted on permanent glycerin slides. Nematodes were identified to the lowest possible taxonomic level (genus/species) using the NEMYS database (Steyaert et al., 2005) and literature therein and the Linnaean Society keys to Nematoda (Platt & Warwick, 1983, 1988; Warwick et al., 1998). Wieser’s (1953) classification was used to distinguish four trophic groups: selective (1A) and non-selective (1B) deposit feeders, epistrate feeders (2A), and predators/omnivores (2B).

Diversity measures were calculated based on the nematode abundance data in the integrated sediment column (0–15 cm). As a fixed number of individuals was identified, several different diversity measures were calculated in order to compare species richness and diversity between treatments. Diversity was expressed by the Margalef’s species richness (d), Pielou’s evenness (J′), Shannon–Wiener diversity index (H′) and the rarefaction index ES(x) (expected number of species). One knot of 30 was used [ES(30)] (i.e. the lowest number of nematodes recorded in one sample) to allow comparisons between different treatments. In order to compare diversity among treatments, the number of nematodes used for calculating the diversity indices was standardized to 30% of the total number of nematodes in the microcosm (i.e. the minimum percentage of nematodes sorted for any core). This was done by randomly selecting 30% of individuals from all nematodes recorded in all depth layers in a given microcosm.

Statistical analysis

The effects of treatment on the total densities of meiofauna and selected major taxa in the integrated sediment column were studied by one-way PERMANOVA with the factor Treatment (TR) with four fixed levels (Control, Macoma, Mya, Cerastoderma). The same procedure but with five fixed levels (Control, Field, Macoma, Mya, Cerastoderma) was used to test for differences in nematode densities between the experimental treatments and field community. Univariate data analyses were performed on Euclidean distance based resemblance matrices using unrestricted permutations of raw data. Multivariate analyses were performed using a three-way crossed PERMANOVA. The experimental design consisted of three factors: Treatment (TR; fixed, with four levels), sediment Depth (DE; fixed, with seven levels) and Replicate (RE; random, with three levels) nested within Treatment (TR). The application of PERMANOVA (permutational ANOVAs that can be used as univariate ANOVAs with P-values obtained with permutation) and nesting the replicates within treatment is a way to deal with the lack of independency of data since the different slices originate from the same core (Braeckman et al., 2011a, b; Urban-Malinga et al., 2014). The analyses of effects of treatment on the vertical profiles of total meiofauna and nematode densities, nematode community structure (composition and densities of nematode species), and organic carbon, TN, and pore water nutrient contents were performed on both raw (to take differences in species/genus abundance into account) and standardized (to account for differences in total nematode abundances among treatments) untransformed, square- and fourth-root transformed nematode genera abundance data (to discriminate between the effects of bivalves on more common and rare genera). Euclidean distance and Bray–Curtis-based resemblance matrices were used for abiotic and biological data, respectively. The Bray–Curtis similarity measure is undefined for two empty samples (sediment depth slices with no meiofauna), and therefore, the zero-adjusted Bray–Curtis resemblance matrix was used, in which a dummy species is added to all samples in the original abundance matrix. Abiotic data (organic carbon, TN and pore water nutrient contents) were normalized prior to analysis. PERMANOVA analyses were conducted using 9999 permutations of residuals under a reduced model. Significant interaction effects were further investigated using a posteriori pairwise comparisons of factor TR within levels of TR × DE. Pairwise tests were based on P values calculated using the 9999 Monte-Carlo permutations procedure [i.e. P(MC)]. Homogeneity of multivariate dispersion across groups was tested prior to a posteriori pairwise comparisons by means of PERMDISP. Distances of group members to the group centroids were tested by permutation within RE (TR) groups (averaged Depths) and in TR × DE groups (averaged Replicates).

Non-metric multi-dimensional scaling (MDS) ordinations were used to visualize the similarities between the treatments and replicates. Analysis of similarity percentages (SIMPERs) was performed to determine the contribution of individual species to the average Bray–Curtis dissimilarity between treatments. SIMPER and PERMANOVA analyses were carried out using the software package PERMANOVA+ for PRIMER (Anderson et al., 2008).

Results

Survival of macrofauna

All but three macrofauna specimens initially added to the microcosms were alive at the end of the experiment. The exceptions were two Cerastoderma individuals and one Macoma specimen, which were found dead at the sediment surface after 2 weeks of incubation and replaced.

In control cores, an oxidized surface zone of 10–15 mm was evident, below which the sediment was grey or greyish-black. In the Cerastoderma treatments, the depth range of the oxidized zone was not visibly different, whilst in the Macoma treatments, it was always extended to a depth of 4–5 cm. Microcosms with Mya either had no clearly visible changes to the oxidized zone or it extended much deeper, up to 10 cm depth. Slicing sediment cores from the microcosms at the end of the experiment demonstrated that Macoma was concentrated mainly in the upper 6 cm, whilst Mya was noted in the upper 10 cm and Cerastoderma in the upper 2 cm.

Sediment characteristics

Average TOC values were enhanced in the upper 4 cm in the Mya treatment, whilst TN contents were lower in the top sediment layer in the Macoma treatment compared to the control and other treatments (Fig. 1). Pore water ammonium concentrations in surface sediments were lower in Macoma and Cerastoderma microcosms compared to the control, whilst in the Mya treatment ammonium content peaked in the topmost sediment and was reduced in deeper sediment layers (Fig. 2). PERMANOVA showed, however, no effect of treatment on vertical profiles of pore water nutrients concentrations, organic carbon or TN contents (Table 1).

Fig. 1
figure 1

Vertical profiles of organic carbon and total nitrogen content in each microcosm treatment at the end of the experiment (mean ± SE; control: empty dots, treatments with bivalves: filled dots)

Fig. 2
figure 2

Vertical profiles of pore water ammonium and phosphates concentrations in each microcosm treatment at the end of the experiment (three replicates integrated; control: empty dots, treatments with bivalves: filled dots)

Table 1 Results of PERMANOVA analyses for differences in vertical profiles of pore water nutrients and organic carbon and total nitrogen contents, total meiobenthic and nematode abundances, multivariate nematode community structure and selected nematode abundances among treatments and sediment depths at the end of the experiment

Effects on meiofauna

Total abundances

Average total abundances of meiofauna (integrated over the sediment column) recorded in the Mya and Cerastoderma microcosms at the end of the experiment were higher than those in the sediment devoid of macrofauna [means (±SD): 658 (±360) and 587 (±156) vs. 412 (±47) ind. 10 cm−2, respectively], and lowest in the Macoma treatment [179 (±112) ind. 10 cm−2] (Fig. 3). These differences were not, however, statistically significant (PERMANOVA, Table 3). Nematodes numerically dominated all treatments, ranging from 58% of the total meiofaunal abundance with Cerastoderma to 88% with Macoma. In addition, juvenile bivalves and rotifers were numerically important groups, constituting 11–20% of the total meiobenthic density in treatments with Cerastoderma and Mya, and 0–7% in the Macoma treatment and defaunated sediment, respectively (Fig. 3). There were significant differences in their abundance across treatments (Table 3). In contrast, integrated nematode abundances at the end of the experiment [means (±SD): 158 (±88)–462 (±70) ind. 10 cm−2] were not statistically different across treatments (PERMANOVA, Table 3).

Fig. 3
figure 3

Total meiofauna numbers (mean ± SE) and numerical abundance of the dominant major taxa in each microcosm treatment at the end of the experiment

Nematodes in the experiment versus in the field

Total nematode densities recorded in the microcosms at the end of the experiment were not significantly different from those recorded in the field on the day of sediment sampling [PERMANOVA: MS = 44,252, Pseudo-F = 1.53, P(perm) = 0.26]. However, the nematode assemblages differed significantly [PERMANOVA: MS = 1257, Pseudo-F = 2.4, P(perm) = 0.003, respectively]. Pairwise tests showed that the field community varied significantly from those recorded in experimental treatments [P(MC) < 0.05] except for the Macoma treatment [P(MC) > 0.05]. SIMPER analysis of the presence/absence data revealed that differences in the nematode composition were mainly due to the following genera: Daptonema, Tripyloides, Viscosia, Eleutherolaimus and Calomicrolaimus (data not shown). Viscosia and Eleutherolaimus were present in the field but not recorded in the microcosms, whilst the densities of Tripyloides, Daptonema and Calomicrolaimus were strongly reduced under experimental conditions (Table 2). Integrated nematode abundances at the end of the experiment [means (±SD): 158 (±88)–462 (±70) ind. 10 cm−2] were not statistically different compared to the field community [447 (±65) ind. 10 cm−2; PERMANOVA: MS = 44,252, Pseudo-F = 1.53, P(perm) = 0.261].

Table 2 Average (mean ± SE) abundances (ind. 10 cm−2) of nematode species in the natural environment and in each treatment at the end of the experiment

Nematode community structure

In total, 26 nematode genera including 5 multispecies genera and 21 species were recorded (Table 2). Adoncholaimus thalassophygas (de Man, 1876) followed by Desmolaimus cf. zeelandicus de Man, 1880 dominated all treatments together constituting from 36% of abundance in the Macoma treatment to 58% in the control. Halomonhystera disjuncta (Bastian, 1865) Andrassy, 2006 was co-dominant, representing 3–12% of abundance.

The effect of treatment on total nematode community structure (integrated over depth layers) was significant when standardized square- and fourth-root nematode abundance data were analysed (Table 3), but a posteriori tests showed that no one pair of treatments was responsible for the difference among treatments [P(perm) > 0.05]. Analysis of the MDS plot (Fig. 4) showed that the community in the Macoma microcosms was separated from those in other treatments. There was no effect of treatment on nematode community structure when raw data were analysed.

Table 3 Results of PERMANOVA analyses for differences in total abundances of meiofauna and selected major taxa, nematode community structure and diversity indices in the integrated sediment column among treatments at the end of the experiment
Fig. 4
figure 4

Non-metric multidimensional scaling (MDS) ordinations of treatment similarity based on standardized square-root transformed nematode abundance data in the integrated sediment column at the end of the experiment

Diversity indices [d, J′, H′, ES(30)] were not statistically different across treatments [P(perm) < 0.05] (Table 3).

Vertical distribution

The majority of meiobenthos in the control and in the microcosms with Mya and Cerastoderma were recorded at the sediment surface (87–96%), whilst in the Macoma treatment, 46% of meiobenthos penetrated deeper sediment layers (Fig. 5). Vertical profiles of both total meiobenthic and nematode abundances were significantly affected by the treatment [P(perm) < 0.05, Table 1]. A posteriori tests showed that meiobenthic and nematode abundances were significantly lower in the top sediment layer (0–1 cm) in the Macoma treatment compared to all other treatments [P(MC) < 0.05; data not shown].

Fig. 5
figure 5

Vertical distribution of meiofauna and major meiofaunal taxa (mean ± SE) in each treatment at the end of the experiment

No effect of treatment on the vertical profiles of nematode community structure was identified both when standardized and raw nematode abundance data were analysed (PERMANOVA, P > 0.05; data not shown).

The vertical distribution of these dominant nematodes (Fig. 6) (A. thalassophygas, D. cf. zeelandicus, H. disjuncta) varied significantly across treatments [P(perm) < 0.05; Table 1]. Pairwise tests showed that their abundances in the top sediment layer (0–1 cm) in the Macoma treatment (29.0, 3.7, 8.7 ind. 10 cm−2, respectively) were significantly lower than in other treatments (95–180, 15–50, 9–55 ind. 10 cm−2, respectively). Also other nematode species were one order of magnitude less abundant in the upper sediment in Macoma treatment than in other microcosms. Deeper sediment layers in Mya and Macoma microcosms were penetrated largely by Ascolaimus elongatus (Bütschli, 1874), D. cf. zeelandicus and A. thalassophygas.

Fig. 6
figure 6

Vertical distribution of nematodes (mean ± SE) and numerical abundance of the dominant nematode species in each treatment at the end of the experiment

Discussion

Our experiment demonstrated that the presence and identity of bivalves did affect vertical distribution of meiobenthos and nematodes, and influenced the composition of meiobenthic and nematode community at the sediment surface.

Methods

The effect of sediment sieving, homogenization and stabilization prior to the addition of macrofauna on meiobenthic community in this experiment is discussed in detail by Urban-Malinga et al. (2014). It was concluded that a vertical meiofauna distribution was re-established during the stabilization period, thus achieving a new equilibrium prior to the addition of macrofauna. Nematode densities recorded in the microcosms at the end of the experiment were not significantly different from those recorded in the field on the day of sediment sampling, suggesting no mortality of nematodes under the experimental conditions. Nematode community structure and composition in the treatments were also similar to the those recorded in the natural environment except in the Macoma treatment in which they were significantly different from the field community. With the exception of Viscosia viscosa (Bastian, 1865) de Man, 1890 and Eleutherolaimus spp. which were in low abundance at the study site (8 and 1.7 ind. 10 cm−2, respectively) and did not survive under experimental conditions, all other species recorded at the study site were also found in the microcosms. These observations suggest that our experimental procedure had little effect on nematode community composition, and therefore, recorded differences can be attributed to the presence of the three bivalve species.

Effect on meiobenthos and nematodes

The addition of the three bivalve species to the experimental microcosms had no significant effect on total meiofaunal or nematode density. There were, however, significantly reduced numbers of two major taxa, juvenile bivalves and rotifers, in the Macoma treatment compared to treatments with other bivalves. Only two specimens of juvenile bivalves were recorded in microcosms with Macoma (vs. 50–106 ind. 10 cm−2 in microcosms inhabited by other bivalves and 14–32 ind. 10 cm−2 in defaunated sediment), whilst the numbers of rotifers were similar to those recorded in defaunated sediment (on average 6 and 3 ind. 10 cm−2, respectively vs. 70–117 ind. 10 cm−2 in microcosms with other bivalves). Juvenile bivalve numbers were probably elevated in Cerastoderma and Mya microcosms owing to the presence of biofilms: faeces and pseudofaeces deposited at the sediment surface by these two bivalves presumably promoted the microbial community which in turn supported the microbial-feeding rotifers. In contrast, the mechanical disturbance of the sediment surface by Macoma activity and surface deposit-feeding, and its release of faeces to the subsurface sediment, were probably responsible for reduced numbers of these taxa in the Macoma treatment.

The observed reduction, although not significant, of the numbers of all other meiobenthic taxa in all Macoma microcosms compared to other treatments must be also noted. Interestingly, this decrease is similar for all the dominant major taxa (i.e. an average drop of 54–66% for Nematoda, 47–72% for Turbellaria and 47–68% for Harpacticoida compared to other treatments). This may indicate meiofaunal mortality either due to competition for food with Macoma or physical disturbance by Macoma feeding activity and general movement. Reduced meiofaunal numbers in response to the activity of M. balthica were reported by Olafsson et al. (1993, 2005), whilst Reise (1981, 1983) observed stimulation of turbellarians in the presence of Macoma and attributed this to organic enrichment of the subsurface sediment (2–4 cm) where the Macoma exhalant siphon probably terminated. In our experiment, organic carbon contents were non-significantly enhanced in the subsurface sediment in Macoma treatment. They were also higher in the upper 4 cm in the Mya treatment than in the defaunated sediment and other bivalve treatments. Hansen et al. (1996) suggested that the large “burrow” openings of M. arenaria might trap labile organic matter when the siphons close or retract, whilst repeated extensions and withdrawals of the long siphons may cause transport of the overlying water with associated organic suspension, faeces and pseudofaeces into the burrow, leading to enrichment of the surrounding sediment. All changes in vertical sediment characteristics observed in our experiment were, however, not significant, suggesting too short incubation time and/or too small number of replicates used to study changes in the sediment parameters (e.g. pore water samples were pooled together to obtain sufficient amount of pore water for analysis). Due to the same reasons, the overall effects of bivalves on meiobenthic and nematode communities, especially integrated over the sediment layers, were weak. It must be highlighted, however, that both time period and replicates number were sufficient to observe significant changes in the vertical profiles of nematode occurrence, especially in the top sediment layer. Nematode occurrence in Cerastoderma microcosms was largely limited to the sediment surface whilst in the Macoma and Mya treatments nematodes penetrated deeper sediment layers. The lack of response of meiobenthos to the presence of Cerastoderma was also reported in earlier studies (Reise, 1983; Kennedy, 1993) and can be related to the near-surface activity of this bivalve, having limited effect on subsurface sediment characteristics. Van Colen et al. (2012) observed a substantial effect of Cerastoderma on benthic processes, but these authors studied muddy sediments (vs. medium sand in our study) and larger North Sea specimens which potentially have a stronger effect on sediment characteristics.

In Mya microcosms, the majority of nematodes were also concentrated in the top sediment layer but, on average, 13% of the community penetrated to deeper sediment layers. A similar habitat extension was observed in Macoma microcosms, where nematodes were recorded to a depth of 10–15 cm and only 54% of the total meiofaunal abundance was limited to the top sediment layer. Similarly, Reise (1981, 1983) observed penetration of turbellarians to deeper layers of sandy sediments in the presence of Macoma. By contrast, in muddy sediment Olafsson et al. (1993) found no effect of Macoma on meiofauna distribution, but in such habitats, the typically severe chemical gradients are likely to be more important than bivalve activity in structuring interstitial communities.

We hypothesize that in our experiment, Macoma and Mya have also created favourable conditions for non-selective (A. elongatus, D. cf. zeelandicus) and surface deposit-feeding [Paracanthonchus spp., Paracyatholaimus proximus (Bütschli, 1874) Micoletzky, 1924)] nematodes, which penetrated deeper sediment layers in these microcosms compared to the defaunated sediment and the Cerastoderma treatment (data not shown).

Differences in the vertical distribution of nematode abundance and nematode species in the Macoma microcosms did not result, however, in significant changes in the vertical structure of the community. Results of pairwise tests showed that abundance changes in the top sediment layer (0–1 cm) were in fact responsible for differences in meiofaunal and nematode assemblages between treatments. The reduced numerical abundance of A. thalassophygas and D. cf. zeelandicus observed in the top sediment layer in Macoma microcosms, compared to all other treatments (including the field community), may be partly due to competition for food between these dominant nematodes and the bivalve and/or due to changes in food availability in response to bivalve mechanical disturbance (Braeckman et al., 2011a; Van Colen et al., 2012). Also, the abundance of the bacterial feeder, H. disjuncta, was significantly reduced in the presence of Macoma in comparison to other treatments probably due to sediment surface disturbance and the absence of a biofilm. These observed differences suggest that the release of pseudofaeces and faeces to the sediment surface by Cerastoderma and Mya and the associated creation of biofilms that facilitate the development of the microbial community may be a significant factor influencing the structure of the meiobenthic community, both at the meiobenthic higher taxa (juvenile bivalves and rotifers) and nematode species levels.

The densities of Macoma used in our experiment were higher than average natural densities of this bivalve, but the biomass corresponded to average macrobenthic biomass in the field. The fact that the effect of Macoma on the occurrence of rotifers, bivalve juveniles, and on nematodes (specifically the reduced abundance of A. thalassophygas, D. cf. zeelandicus and H. disjuncta) observed here is similar to that observed when another active surface deposit feeder, the polychaete Hediste diversicolor (O. F. Müller, 1776), is studied under the same experimental conditions (Urban-Malinga et al., 2014) must be highlighted. This observation confirms the prominent role of macrobenthic functional traits in structuring meiobenthic communities. Since the bivalves described here co-occur in the natural environment, it would be interesting to investigate their interactions. There was no indication of any effects between sessile M. arenaria and C. glaucum under experimental condition (Urban-Malinga et al., 2014) but M. balthica living actively closer to the sediment surface than M. arenaria and feeding at the sediment surface is likely to interact with surface-dwelling C. glaucum.