Evaluating environmental drivers of spatial variability in free-living nematode assemblages along the Portuguese margin

Understanding processes responsible for shaping biodiversity patterns on continental margins is an important requirement for comprehending anthropogenic impacts in these environments and further management of biodiversity. Continental margins perform crucial functions linked to key ecological processes which are mainly structured by surface primary productivity and particulate organic matter flux to the seafloor, but also by heterogeneity in seafloor characteristics. However, to what extent these processes control local and regional biodiversity remains unclear. In this study, two isobathic parallel transects located at the shelf break (300–400 m) and upper slope (1000 m) of the western Iberian margin were used to test how food input and sediment heterogeneity affect nematode diversity independently from the spatial factors geographical distance and water depth. We also examined the potential role of connectedness between both depth transects through molecular phylogenetic analyses. Regional generic diversity and turnover were investigated at three levels: within a station, between stations from the same depth transect, and between transects. High variability in food availability and high sediment heterogeneity at the shelf-break transect were directly linked to high diversity within stations and higher variation in community structure across stations compared to the upper slope transect. Contrastingly, environmental factors (food availability and sediment) did not vary significantly between stations located at the upper slope, and this lack of differences were also reflected in a low community turnover between these deeper stations. Finally, differences in nematode communities between both transects were more pronounced than differences within each of the isobathic transects, but these changes were paralleled by the previously mentioned environmental changes. These results suggest that changes in community structure are mainly dictated by environmental factors rather than spatial differences at the western Iberian margin. Furthermore, phylogenetic relationships revealed no evidence for depth-endemic lineages, indicating regular species interchanges across different depths.


Introduction
The link between biodiversity (i.e. species diversity) and ecological processes (e.g. carbon flow, surface productivity) has created a heightened interest in ecological research after large-scale human impacts were deemed responsible for declining species numbers and alterations of ecosystem properties (Loreau et al., 2001). Stretching between the coastal zone and the abyssal plains of the deep sea, continental margins (100-4000 m) encompass the largest habitat diversity in the marine environment (Levin and Dayton, 2009; Ramirez-L. Lins et al.: Evaluating environmental drivers of spatial variability in free-living nematode assemblages are responsible for 90 % of the new biological productivity in oceans and seas, providing valuable food and energy resources for the marine fauna (Salgueiro et al., 2014).
It is generally accepted that principal biological oceanographic processes, such as carbon burial and nutrient cycling, remain concentrated within continental margins (Levin and Dayton, 2009). However, the biodiversity of continental margins is under severe threat by increasing commercial exploitation, ranging from fisheries to gas, oil, and mineral extraction (Levin and Dayton, 2009;Puig et al., 2012). The direct impact of these unabated commercial activities on the benthic environment and populations varies greatly from pervasive sediment erosion, transportation, and deposition to the large-scale alteration of community composition (Puig et al., 2012). Therefore, continental margins comprise key locations to study the effects of environmental alterations on benthic biodiversity.
Understanding the processes that shape biodiversity patterns on continental margins is an important prerequisite for comprehending and managing anthropogenic impacts in these environments. Sea-surface processes have an important effect on the benthic fauna because part of the primary production is exported from overlying waters to the deep-sea floor, mostly in the form of phytodetritus, where it serves as food source to benthic communities (Billett et al., 1983;Lins et al., 2015;Serpetti et al., 2013;Wei et al., 2010). Particulate organic carbon input in the deep sea has been regarded as one of the main factors shaping benthic community structure and functioning (Rex, 1981). Phytodetritus creates patchiness, enhancing habitat heterogeneity, and consequently promotes species coexistence (Cardinale et al., 2000). In addition, water depth indirectly plays a role in structuring benthic communities, since organic matter flux is negatively related to depth, and deeper regions will consistently receive less input of labile organic matter compared to shallower regions Garcia and Thomsen, 2008;Lutz et al., 2007;Ramalho et al., 2014). As a consequence of this decline in food availability, decreases in abundance and biomass associated with an increase in depth on the continental slopes have been observed for most benthic size classes (mega-, macro-, and meiofauna) (Flach et al., 2002;Muthumbi et al., 2011;Rex et al., 2005;Rowe et al., 2008;Thiel, 1978).
Food availability, as well as biological factors (predation, competition, dispersal), are assumed to drive small-scale (1-10 m) patterns of benthic communities (Gage, 1997) promoting local diversity (Levin et al., 2001). Beta diversity over large spatial scales (100-1000 m) within continental margins have been observed both along bathymetric gradients as well as between different sites at similar depths. This indicates that beta diversity is not singularly depth dependent (Danovaro et al., 2013;Easton and Thistle, 2016;Havermans et al., 2013;Leduc et al., 2012a). Physical factors, including near-bottom currents, sediment grain-size heterogeneity, boundary constraints, human activities, and topography, are also considered of particular importance for beta diversity (Levin et al., 2001). They shape biodiversity as they may reduce the effect of a dominant species through the redistribution of resources among inferior and superior competitors (Stachowicz et al., 2007), in this way increasing species diversity.
Moreover, population dynamics and dispersal (Derycke et al., 2013;Gage, 1997;Rex et al., 2005) have been shown to affect the structuring of benthic fauna at different spatial scales. Most benthic species have restricted active dispersal potential, although passive dispersal may be facilitated through ocean currents (Etter and Bower, 2015;Gallucci et al., 2008;Lins et al., 2014;Ullberg and Olafsson, 2003). The lack of a pelagic larval stage in free-living nematodes, the focus group of this study, could therefore be viewed as a disadvantage to dispersal. Nevertheless, this abundant and omnipresent group of benthic metazoans is found at all depths and in all deep-sea habitats (Giere, 2009;Vincx et al., 1994). Nematodes, which belong to the meiofauna (< 1 mm), exhibit high species richness and are one of the few taxa in which true cosmopolitan species may exist (Bik et al., 2010;Zeppilli et al., 2011). Some species are able to actively swim, following chemical cues, but more importantly, nematodes may be passively transported via water currents following resuspension from disturbance events (Jensen, 1981;Schratzberger et al., 2004). Molecular studies have indicated that different nematode taxa in diverse habitats exhibit population connectivity across a wide geographical range, with some species showing subtle but significant genetic structuring at a small spatial scale and other species exhibiting no differentiation along large distances (> 500 km). These findings confirm a high dispersal potential and low endemism for at least some species (Derycke et al., 2013(Derycke et al., , 2005. Nematodes therefore hold ideal life history traits when seeking to understand dispersal, coexistence, and benthic-pelagic coupling in the deep sea. Depth-related factors are thought to inhibit across-depth gene flow and thus to promote speciation in some taxa. This depth-range limitation provides another explanation for why the bathyal holds such a high biodiversity (Rex and Etter, 2010). However, while empirical data for macrofaunal molluscs, crustaceans, as well as octocorals have been found supporting this depth-differentiation hypothesis (France and Kocher, 1996;Jennings et al., 2013;Quattrini et al., 2015), it may not apply to nematodes, for which repeated and regular interchanges between depths have been observed (Bik et al., 2010).
Most previous research on nematode diversity has concentrated either on bathymetric differences (Danovaro et al., 2013;Muthumbi et al., 2011) or on geographical transects and macro-habitat heterogeneity (Baldrighi et al., 2014;Lambshead et al., 2000;Van Gaever et al., 2009). No studies so far have combined bathymetric with geographic analyses along an isobathic transect at a regional scale, which is crucial to understand patterns of biodiversity. Furthermore, the transition from the shelf to the slope remains a largely understudied area (Muthumbi et al., 2011;Vanreusel et al., 1992), while major environmental shifts are observed here. Through the analysis of two isobathic transects of about 20 km length, and separated by 30 km and 600 m water depth, we tested the effect of environmental (food availability and sediment heterogeneity) and spatial (depth and geographical distance) variables on nematode diversity. In this way, drivers for turnover in nematode taxonomic composition were analysed at three spatial scales: within stations (i.e. beta diversity between cores from the same station), beta 1, between stations from the same depth (i.e. beta diversity between stations belonging to the same depth transect), beta 2, and between two depth transects, beta 3. As food input and sediment heterogeneity are expected to vary more with depth than along regional isobathic transects, we expected a higher turnover in community composition between depths than within isobathic transects. To evaluate possible depth-mediated differentiation, genus turnover and phylogenetic relationships through DNA sequence clustering between the two isobathic transects were investigated based on 18S rDNA sequence data of selected nematode taxa. The following hypotheses were tested: (H1) higher patchiness of food resources deposited at the seafloor results in a higher diversity at small spatial scales (i.e. beta 1); (H2) increased sediment heterogeneity results in a higher beta diversity; (H3) beta diversity between different bathymetric transects is higher than beta diversity across similar depths; (H4) clades/taxa are shared between shelf break and slope areas.
2 Material and methods

Sampling and study area
The western Iberian margin (WIM) is characterized by a narrow shelf and steep slope (Garcia and Thomsen, 2008;Nolasco et al., 2013;Relvas et al., 2007). Primary production in this area increases in May-June and constitutes a significant proportion of the yearly production, reaching values higher than 90 g C m 2 y −1 (Salgueiro et al., 2014). The WIM exhibits seasonal upwelling with filaments that can penetrate more than 200 km into the open ocean, influencing not only vertical transport but also horizontal particle transport from near shore towards the open ocean (Crespo et al., 2011;Figueiras et al., 2002;Relvas et al., 2007;Salgueiro et al., 2014Salgueiro et al., , 2010. The high particle transport observed at the WIM occurs mainly due to the great bottom dynamics in the area. This region possesses an equatorward current flow generated by thermohaline structures of water masses and wind forcing, eddy interactions with the alongshore circulation, and buoyant plumes (Relvas et al., 2007). These features, together with shelf and coastal currents, upwelling filaments, and fronts, impact the subsurface circulation, internal waves, and consequently the transport of sinking particulate organic  , 1997;Relvas et al., 2007). During the RV Belgica B2013/17 (10-18 June 2013) and B2014/15 (2-10 June 2014) cruises to the WIM, sediment samples for nematode and environmental analyses were taken at the slope off the southwest coast of Portugal (Fig. 1). The study area comprised two main transects roughly parallel to the isobaths. The first transect was 23 km long, situated 294-445 m deep (further referred to as shallow transect), just beyond the shelf break; the second transect was located at the upper slope, 19 km long, and at a water depth of 900-1006 m (named deep transect). The "shallow" transect included six stations while the "deep" transect comprised four stations (Table 1). Sampling was performed using a Multicorer (MUC) equipped with four Plexiglas tubes yielding samples with a virtually undisturbed sediment surface (inner core diameter 9.8 cm).

Sediment analyses
Three to four replicated samples for granulometric and geochemical analyses (1 g of sediment) from the first sediment layer (0-1 cm) were frozen at −80 • C. Grain-size distribution was measured with a Malvern Mastersizer 2000 (0.02-2000 µm size range) and divided into five categories, from silt clay to coarse sand fractions. Sediment particle-size diversity (SED) was calculated as a measure for sediment heterogeneity from the percent dry weight of the five size classes mentioned above using the Shannon-Wiener diversity index (Etter and Grassle, 1992;Leduc et al., 2012b). Total sedimentary organic carbon (% TOC) and nitrogen (% TN) were determined with a Carlo Erba elemental analyser on freeze-dried and homogenized samples after acidification with 1 % HCl to eliminate carbonates. Total organic matter (% TOM) content was determined after combustion of the sediment samples at 550 • C.
Chlorophyll a (Chl a), chlorophyll degradation products, and carotenoids in the sediment were measured with a Gibson fluorescence detector (Wright and Jeffrey, 1997) after lyophilization, homogenization, and extraction in 90 % acetone, and separation of the samples via reverse-phase HPLC (high-performance liquid chromatography). Chloroplastic pigment equivalents (CPE: Chl a + phaeopigments) were used as a proxy for surface-derived primary productivity at the seafloor.

Nematode sample processing for community analyses
At each station, three to four replicate samples of the 0-1 cm layer were used for nematode analysis. Samples were fixed on board with seawater buffered 4 % formalin. Sediment was washed over 1000 and 32 µm sieves. The fraction retained on the 32 µm sieve was centrifuged three times using LUDOX HS40 Dupont (specific gravity 1.19) as flotation medium and then stained with Rose Bengal. In each sample, 140 nematode individuals (whenever enough present) were randomly picked out and gradually transferred to glycerine (De Grisse, 1969), mounted on glass slides and identified to genus level using relevant literature (Vanaverbeke et al., 2015;Warwick et al., 1998). Functional diversity (based on the relative abundance of each trophic type) of nematodes was calculated using individuals trophic levels according to Wieser (1953): selective deposit feeders (1A), non-selective deposit feeders (1B), epistratum feeders (2A), and predators (2B), complementing the 2B group with the notion of "scavengers" (Jensen, 1987). Trophic diversity (TD) was calculated using the index proposed by Heip et al. (1985): where q i is the relative abundance of trophic type i. Taxonomic diversity was measured using Shannon-Wiener diversity (H ), expected nematode genus richness (EG (80)), and Pielou's evenness (J ).

Data analysis
Trends in environmental variables (% TOC, % TN, % TOM, Chl a, CPE, carotenes, depth, sediment grain size, and SED) and univariate nematode variables (H , J , EG (80), and TD) were investigated by means of Spearman rank correlations and Draftsman plots (Anderson et al., 2007) in R (R Core Team, 2013).
The untransformed nematode community data based on relative abundance of genera were analysed based on Bray-Curtis similarities (and Euclidean distances for the univariate data) by means of non-parametric multivariate ANOVA (PERMANOVA; Anderson et al., 2007) to assess differences between "deep" (slope) and "shallow" (shelf-break) areas (two-factor nested design). The two-factor model included "depth" as a fixed factor and "station" as a fixed factor nested in "depth". Due to the use of an unbalanced design, the type I of sum of squares was chosen for the PERMANOVA analysis to make sure all possible rearrangements of samples are equally likely (Anderson et al., 2007). Subsequent pairwise pseudo t tests were performed between all pairs of levels to determine where significant differences between each combination were found. Additionally, PERMDISP routines were used to test for homogeneity of multivariate dispersions between stations. PERMDISP results were not significant, indicating location differences through equally dispersed distances to centroids. SIMPER routines were executed based on Bray-Curtis similarity, with a cut-off of 90 % for low contributions. Dissimilarities within and between stations were compared with distances between geographical areas (km) and between depth differences (m).
Diversity indices and nematode community structure turnover were investigated at three hierarchical levels: within stations, beta 1 (i.e. beta diversity between cores from the same station), between stations from the same depth, beta 2 (i.e. beta diversity between stations belonging to the same depth transect), and between two depth transects, beta 3. Turnover between stations from the same transect and between transects were statistically tested using the two-factor nested PERMANOVA design described above. Dissimilarities within and between stations and between transects were used to compare variability within each spatial scale.
The multivariate environmental data were first normalized (subtracted mean divided by standard deviation) and resemblance matrices were calculated based on Euclidean distances. Subsequently, PERMANOVA tests were performed using the same design as described for the multivariate nematode community data. PERMDISP analyses were conducted on the multivariate environmental data to investigate variability in the environmental conditions between the two transects. Subsequently, SIMPER routines were executed to identify which factors were responsible for between-transect differences. DistLM (distance-based linear model) routines were performed to analyse and model the relationship between nematode genus community and environmental variables with correlations lower than 0.9 (Chl a, carotenes, CPE, % TN, silt clay, very fine sand, medium sand, and coarse sand). Highly correlated variables (% TOC, % TOM, fine sand, and depth) were first transformed to cosine and if high correlations persisted they were excluded from the DistLM analysis. The DistLM model was built using a step-wise selection procedure and adjusted R 2 as a selection criterion. Euclidean distance was used as a resemblance measure for DistLM procedures and the results were displayed in dbRDA (distance-based redundancy analysis) plots.

Molecular phylogenetic analyses of nematodes
One sample (replicate) from each of the "shallow" stations S4 and S2 and one from the "deep" station D4 were preserved in DESS (Yoder et al., 2006) and used for molecular analyses. The first centimetre (0-1 cm) of each core was washed with LUDOX HS40 Dupont, following the same protocol as for the community analysis (see above). One hundred nematodes were randomly picked out per sample under a stereomicroscope (50× magnification). Each individual was rinsed in sterile water, transferred to a microscope slide containing sterile water, and digitally photographed as morphological reference with a compound microscope Leica DMR and Leica LAS 3.3 imaging software. DNA extraction followed Derycke et al. (2005) using the entire specimens.
Successful PCR reactions were identified using agarose gels stained with ethidium bromide and were sequenced with both forward and reverse primers by Macrogen Europe (the Netherlands) with the fluorescent dye terminator Sanger sequencing method. The resulting reads were assembled using Mega 6.0 (Tamura et al., 2013). Sequences were checked for contamination using the BLAST algorithm on GenBank (Benson et al., 2008). The sequences that showed contamination or were of low quality (high amount of ambiguous nucleotides) were removed from the alignment. Nematode contig sequences (consensus of forward and reverse sequences) generated during this study were aligned using the MAFFT multiple sequence alignment algorithm (Katoh et al., 2009) as implemented in Geneious 9.0 (Kearse et al., 2012) at default settings (the alignment algorithm was automatically determined; scoring matrix was 200PAM/k = 2; gap-opening penalty was 1.53 and the offset value was 0.123).
GenBank sequences for the most representative genera in the samples (all of the nematodes belonging to the class Chromadorea) were included from GenBank (Benson et al., 2008) (whenever available) to compare differences in genetic-/phylogenetic diversity between different depths and locations. Sequences from Meldal et al. (2007) and from Bik et al. (2010) were used to compare genus diversity and diversity within the genus Halalaimus, respectively, between different habitats.
For both datasets, Modeltest 2.1 (Posada and Crandall, 1998) and jModeltest (Posada, 2008) were used to determine that the best suitable model for maximum likelihood (ML) analyses of the nuclear data was according to the Akaike information criterion (AIC) (Akaike, 1981) general time reversible with gamma distributed rate variation among sites and proportion of invariable sites (GTR + G + I).
Reconstruction of 18S relationships was conducted using ML. The analyses were performed by means of randomized axelerated maximum likelihood (RAxML) (Stamatakis, 2006) in raxmlGUI (Silvestro and Michalak, 2012) using the fast likelihood search with 1000 replicates to calculate bootstrap support values.
For the Halalaimus dataset, Bayesian inference was additionally applied in MrBayes (Ronquist and Huelsenbeck, 2003) to supplement topological inferences. Analyses were run for 5 000 000 generations using six MCMC chains. From all runs the first 25 % of sampled trees were discarded as burn-in. Consensus trees were used for illustration here and were ordered and annotated in FigTree and Geneious tree viewer and colorized in Adobe Illustrator. In the supplement tree, line thickness indicates strength of bootstrap support.

Environmental variables
Biogeochemical and granulometric properties of the sediment are shown in Fig. 2. Nested PERMANOVA results for SED revealed significantly lower values (p< 0.05) at the "deep" transect compared to the "shallow" transect. Pairwise comparisons for SED showed significant differences for the pairs [D2, D3] at the "deep" transect and for all the pairs at the "shallow" transect, except for [S1, S4], [S2, S614], and [S2, S613]. Sediment composition at the "deep" stations was mainly composed of silt-clay fractions (81-89 %), while at the "shallow" stations sediment was more heterogeneous. At the shelf-break, fine sand (25-42 %) dominated, except for S7, where medium sand showed a higher proportion (30 %). Nested PERMANOVA results showed significant differences in sediment composition between depth transects and among stations within the same transect (p< 0.05) (Table S1 in the Supplement). Pairwise comparisons between stations showed higher variability in sediment composition for the "shallow" stations, where the pairs of stations [S7, S2], [S2, S614], [S2, S613], and [S614, S613] showed similar sediment characteristics (Table S1). Pairwise comparison for "deep" stations only showed differences between D2 and D3. Within-station comparison showed low variability (< 25 % deviation from the mean values) in silt clay and very fine sand for most stations both shallow and deep (Fig. 2). Fine, medium, and coarse sand variability within each station was higher when compared to silt clay and very fine sand (Fig. 2). Significantly higher values (nested PERMANOVA, p< 0.05) of % TOM (Table S2), % TOC (Table S3), and % TN (Table S4)  No strong variability (< 25 % deviation from the mean values) was observed within stations for these three variables (Fig. 2). Chl a (0-0.17 µg g −1 ), carotenes (0-0.72 µg g −1 ), and CPE (0.01-1.79 µg g −1 ) values were generally low. Chl a showed no significant differences between depth transects (p> 0.05) and only the pairs [S1, S7] and [S1, S614] were significantly different from each other (Table S5). In addition, Chl a showed high variability (> 25 % deviation from the mean values) at the "shallow" stations, especially at S4 (Fig. 2). Carotenes and CPE also possessed high variability (> 25 % deviation from the mean values) at the shelf-break stations and revealed significant differences between depths and among pairs of stations (p< 0.05). For carotenes, the pairs of stations [S1, S4], [S1, S7], and [S1, S614] were significantly different from each other (Table S6), while for CPE only the pair [S1, S613] was significantly different (Table S7). Moreover, carotenes were completely absent at the "deep" stations ( Fig. 2).
TD revealed significantly higher values at the "shallow" stations (p< 0.05) but no significant pairwise difference was observed for the "shallow" or "deep" transect (p> 0.05) (Fig. 4). Nested PERMANOVA results for relative abundance of trophic groups displayed significant differences between the two transects, but not among stations from the same depth (p> 0.05) (Table S12). Average similarity between "deep" and "shallow" areas was 80 %. SIMPER analyses revealed that differences between depths were mainly due to the higher relative abundance of selective deposit feeders (1A) at deeper stations. The shallow stations exhibited higher abundance of epistratum feeders (2A) and predators/scavengers (2B).

Correlation between nematode community structure and environmental variables
The correlation between univariate diversity values (H , J , EG (80), and TD) and environmental variables (% TN, % TOC, Chl a, carotenes, CPE, sediment grain size, and SED) are shown in Table S13. J was not correlated to any environmental factor. Shannon-Wiener diversity (H ) was negatively correlated with silt clay but positively correlated to very fine-medium sand, SED, and total carbon and nitrogen. The EG (80) was negatively correlated with % TN, % TOC, and silt clay and positively correlated with CPE ( Fig. 6), very fine sand, fine sand, medium sand, and SED (Fig. 6). TD was negatively correlated with % TN and silt clay and positively correlated with Chl a, CPE (Fig. 6), very fine sand, fine sand, and SED (Fig. 6). PERMDISP results revealed a higher variability (2.4 ± 0.35 vs. 0.75 ± 0.087, after normalization of the data) in environmental conditions for the "shallow" transect when compared with the "deep" transect. According to the SIMPER results, these differences were mainly derived by higher average values of silt clay and %TN in the "deep" transect, and by the higher coarse sand content in the shallow transect, which together accounted for 45 % of the total dissimilarity. DistLM analyses based on twelve environmental variables explained 33 % of the total nematode diversity. Silt clay accounted for 23 % of the total variation, being responsible for most differences found between "shallow" and "deep" stations (Fig. 7a). The other variables did not contribute significantly to the model and/or added < 5 % in explaining the total variation. When only "shallow" stations were included in the model, the significant environmental variables explained 23 % of the total variation (Fig. 7b). Coarse sand was the main factor accounting for variation between "shallow" stations (13 %). This sediment fraction showed highly fluctuating values (0.02-14.88 %) between stations of the "shallow" areas.

Nematode molecular phylogenetic analyses
From the 300 vouchered nematodes, the success rate of sequencing was only 30 %. For 199 specimens no PCR product was detected or sequences were of low quality. Phylogenetic analyses showed that the 101 sequenced nematodes belong to seven different orders of free-living marine nematodes (Table  S14). The highest genetic diversity was reported for the order Enoplida, with 25 different 18S sequences, followed by the order Plectida (19 different 18S sequences) and Desmodorida (18 different 18S sequences). The ML phylogeny inferred from 18S sequences is shown in Fig. S1 in the Supplement.
In general, the backbone of the Chromadorea phylogeny was poorly supported, leading to several paraphyletic or polyphyletic orders and some families, such as Plectida, Desmodorida, and Oxystominidae. Well supported were the orders Tylenchida (bootstrap support (bs) = 99), Monhysterida (bs = 100), Dorylaimida (bs = 100), Monochida (bs = 100), and Triplonchida (bs = 82). Desmodorida is polyphyletic in our analysis with the family Microlaimidae forming a well-supported clade (bs = 100). The orders Chromadorida and Enoplida represent monophyletic but extremely weakly supported groups (bs = 5 and 41, respectively), while the orders Trefusiida and Triplonchida appeared nested within Enoplida. However, resolving the phylogenetic ties within Chromadorea was not within the scope of this article. What the consensus shows is that the 18S phylogeny supported the broad taxonomic representation of nematodes in the samples and furthermore indicated neither geographic nor depth clustering between "deep" and "shallow" taxa at any level of the tree topology (Fig. S1). This was moreover demonstrated within the best-represented and monophyletic (Fig. S1) genus in the dataset, Halalaimus (15 individuals, complemented with 42 GenBank sequences from different depths and locations globally distributed). Here, the new sequences showed no clustering related to depth or geography was observed but instead they seem randomly scattered between samples from different depths and regions (Fig. 8).

Discussion
Beta diversity within and between stations (beta 1 and 2) from the same transect along the WIM were higher at the shelf break ("shallow" stations) when compared with the upper slope ("deep" stations). The lower beta 1 and beta 2 turnover observed at the upper slope occurred because nematode communities were in general dominated by the same genera (Halalaimus and Acantholaimus), while at the shelf break not only was a higher beta 1 diversity present but also the genus composition clearly differed between stations (beta 2), resulting in a higher turnover. This higher beta 1 and beta 2 diversity and along transect turnover at the shallower stations coincided with a higher amount and patchiness in food supply, as well as a higher sediment heterogeneity within and between stations. However, the largest community differences and therefore highest taxonomic turnover in Figure 5. Dissimilarity values in nematode genus composition among geographical distance (a), sediment composition (b), and chloroplastic pigment equivalents (CPE). Green squares represent dissimilarities between shallow sites, dark blue circles between deep sites, and light blue triangles show dissimilarities between shallow and deep. nematode genera was present between transects (beta 3), supported by the higher dissimilarity values observed compared with within-transect variability. Nevertheless, the high number of shared genera between transects and the intermingled pattern of genetic clustering observed for Halalaimus suggested that depth differences did not restrict the distribution of nematodes.

Sediment heterogeneity in combination with
increased amount and patchiness of food resources contribute to a higher beta diversity In this study, nematode beta 1 and beta 2 diversity was significantly correlated with food availability. Moreover, the higher amount and variability (i.e. dissimilarities within and between stations from the same transect) of "labile" food resources (Chl a, CPE, and carotenes) within the "shallow" stations when compared to the "deep" stations ( Fig. 2), were associated with a higher structural and trophic diversity. The interactions between coastline features and wind forcing already reported for the WIM certainly affect the export flux and, consequently, the resource distribution on the seabed. For this region, the strong surface dynamics may explain the high variability of environmental variables, such as patchiness of food input and habitat heterogeneity, observed here, and consequently the high nematode beta 1 and beta 2 diversity exhibited at the "shallow" stations (Cardinale et al., 2000;Crespo et al., 2011;Tokeshi, 1999). Concurrently, a different pattern was observed for the deeper transect, where the general lower amount and within and between-station variability of food resources were associated with a lower nematode beta 1 and beta 2 diversity. A similar relationship between decreasing availability of "labile" organic matter associated with a decrease in diversity was also observed at similar depths in other slope environments (Danovaro et al., 2013;Leduc et al., 2012a;Netto et al., 2005) and may thus represent a general pattern at continental margins. Nevertheless, the relationship between food and beta diversity in the deep sea is still not well understood and seems to vary according to the different habitats studied. It appears that typical deep-sea communities dominated by similar genus composition can coexist in a food-depleted environment, while seemingly contrasting, biodiversity at both local and regional scale increases with increasing food input. In this sense, when food input reaches a certain threshold, more competitive/opportunistic species are going to proliferate and dominate the community, setting a limit to the diversity of the system Whittaker et al., 2001). This process may possibly explain why diversity in the deep sea is in general higher (at the slope) compared with coastal sediments (Gray, 2002;Moens et al., 2014), despite the higher food input in shallow waters.
However, the lower beta 1 and beta 2 diversity at the slope stations compared to the shelf break is not only explained by the lower and less patchy food input. In this study, beta 1 and beta 2 diversity increased with increasing sediment dissimilarity and SED. Furthermore, SED showed positively correlated values with the expected number of genera and the trophic diversity, exhibiting even stronger correlation results than the one between nematode diversity and CPE (Fig. 6). In addition, DistLM and SIMPER results based on the environmental data revealed a major significant effect of the silt-clay proportion on community beta 3 diversity (Fig. 7).  The high silt-clay contribution, together with the low variation in sediment composition within and between stations of the "deep" transect, was associated with a lower beta 1 and beta 2 diversity in this area (Fig. 5). In addition, the lower TD and higher contribution of deposit feeders in the deep stations were potentially caused by a combination of low food input and comparatively higher sediment stability for this transect, reflected by the finer sediment composition in Figure 8. Consensus of Bayesian inference of phylogeny of the genus Halalaimus based on 18S rDNA sequence fragments generated in this study and from Bik et al. (2010); node support is given as posterior probabilities (PP) and nodes with PP smaller than 0.50 were collapsed. The outgroup was set to Wieseria. The tree shows multiple instances of close relationships between individuals collected at different depth zones. comparison with the "shallow" stations (Fig. 2). The stability of the "deep" transect seems to reduce the generic diversity, and favours the dominance of the genera Acantholaimus and Halalaimus, commonly abundant in deep-sea soft sediments . This decline in TD with depth, positively associated with a decline in SED, was not observed in other studies along the slope, although the SED and TD val-ues observed here were much higher than already reported before for other areas and might reflect different environmental conditions at the WIM (Danovaro et al., 2013;Leduc et al., 2012b;Pape et al., 2013).
In this regard, the homogeneous silty sediments and the persistent low food availability along the "deep" transect might reduce environmental variability and may conse- Figure 9. Scheme showing (a) how beta diversity varied across stations and between bathymetrical transects. The green bar represents the "shallow" transect and the light blue bar the "deep" transect. Grey circles inside the bars represents the stations sampled at each transect. Coloured circles inside grey circles refer to the variability in nematode genus composition within station and across stations. It illustrates the higher densities and turnover (beta 2) found at the "shallow" stations and the lower turnover (beta 1) found at the "deep" stations. (b) The figure reveals the main environmental factors responsible for beta diversity between both depth transects. Upwelling effects, primary production, currents, and decrease of POC (particulate organic carbon) with increasing depth are considered to be correlated with turnover. The fate of organic matter produced at the surface varies with depth, where deeper areas will receive lower labile organic matter (OM) when compared to shallower areas. This figure also illustrates the higher variability in sediment composition at the "shallow" transect when compared to the deeper transect.
quently limit the genus turnover within and between stations there. Contrastingly, the higher within and between-station TD at the "shallow" stations suggests a higher niche differentiation within this transect and the possible transient stage of this habitat. Moreover, higher prominence of opportunistic species and fast colonizers, such as Microlaimus, was observed at the shallower WIM stations (Fig. 3). The high abundance of opportunists has been reported for several bathyal areas in association with disturbance events (Muthumbi et al., 2011;Pape et al., 2013;Raes et al., 2010), pointing to potentially ongoing recolonization processes following the disturbances (Lee et al., 2001). The observed high densities of Microlaimus, which is considered both tolerant to disturbance and an early colonizer (Lee et al., 2001;Moreno et al., 2008;Raes et al., 2010), are in accordance with the assumed hydrodynamic regimes at the shallower stations and possible anthropogenic disturbance effects, for instance from fisheries, which further influence the bottom dynamics there (ICES, 2008;Quaresma et al., 2007;Relvas et al., 2007). Disturbance effects via either bedload movement or erosion and sedimentation of suspended load alters not only particle size but also organic content. Thus, sediment heterogeneity in combination with increased patchiness of food resources deposited at the seafloor observed here were associated with a higher beta 1 and beta 2 diversity and trophic diversity at the shelf break compared with deeper areas. Disturbance has been suggested before as a strong driver for diversity in the deep sea, with disturbances ranging from smallscale bioturbation traces, to intermediate-scale phytodetritus falls or large-scale currents (Levin and Dayton, 2009). In our study, similar disturbances, resulting in higher sediment heterogeneity and food patchiness at the stations located at the shelf break, are used to explain the higher beta 1 and beta 2 diversity compared with the 600 m deeper transect. Actually, the higher beta 1 and beta 2 diversity at the shelf break compared to the upper slope seems to contrast with the regular observed bathymetric diversity gradients for the deep sea, with mid-slope diversity maxima recorded for multiple taxa (Levin and Dayton, 2009;Rex, 1981). However, the higher diversity at genus level for nematodes at the shelf-break stations does not necessarily imply a higher species diversity for this environment, as some of the dominant deeper water genera such as Acantholaimus and Halalaimus are known as highly diverse genera (De Mesel et al., 2006;Muthumbi and Vincx, 1997). Still, generic diversity may better represent the actual functional diversity for nematodes, as different genera are assumed to differ more in function than species within the same genus (Pape et al., 2013).

Beta diversity between different bathymetric transects is higher than beta diversity across similar depths
In general, dissimilarity in nematode genus structure increased with increasing geographical distance and DistLM results showed no significant effect of water depth on the nematode community structure. However, genus turnover between the "shallow" stations was greater than between the "deep" stations. If geographical distance was an important factor, one would also expect a similar high turnover for the "deep" transect as observed for the "shallow" transect as both transects are equally long. Concerning water depth, other studies on various taxa have shown that even small bathymetrical changes can be more important for promoting taxonomic differentiation than large geographical distances within the same depth (Havermans et al., 2013;Quattro et al., 2001). Nevertheless, water depth did not have a significant effect on the differences observed for the nematode community structure in this study and explained only 3 % of the total nematode community variation between transects, appointing environmental differences as the main responsible factors for the higher beta 3 diversity rather than spatial differences. These findings support the general idea that diversity changes can be associated with both large and small-scale features driven by environmental alterations (Fig. 9).
In general, beta diversity in the deep sea appears to be regulated by mechanisms of energy availability, biological interactions, disturbance, and habitat heterogeneity (Levin et al., 2001). Changes in these features at the slope, like the WIM study area, are more pronounced with increasing water depth than with increasing isobathic distances (Rex, 1981), even for short distances. Here, however, a considerable percentage of shared genera still occurred at both transects and some genera exhibited similar abundances for both shelf break and upper slope, such as Tricoma, Daptonema, and Halalaimus (Fig. 3). In this respect, our results do not support the idea of isolation by depth between the shelf break and slope at the WIM area (Bik et al., 2010;Riehl and Kaiser, 2012). To conclude, heterogeneity in sediment and possibly food availability to a lesser extent explain the main depth turnover patterns in this study, while depth and geographical distance are not the main cause for variations in community composition.

Absence of depth-specific clades
Although we observed differences in community structure between the shallow and deep stations, the large proportion of genera shared between the two depth transects may indicate that bathymetrical isolation between the respective populations does not exist (Fig. 3). However, such depth differentiation could occur on the intra-specific and species level. In our study, phylogenetic relationships within and between the genera sampled at both shallow and deep stations revealed shared clades (Fig. S1), potentially representing eurybathic species and thus at least some degree of connectivity across depths, although we can make no conclusions about spatial or timescales. The precise understanding of spatial variability and the processes which drive species diversity and connectivity in the deep sea are still poorly understood (Danovaro et al., 2013;Etter and Bower, 2015). Deep basins are confluent at extensive depths and connected by thermohaline circulation, suggesting they do not represent completely isolated systems (Levin et al., 2001). Processes such as deep-water formation and upwelling, potentially represent means of (passive) across-depth dispersal (Brandt, 1992;Brandt et al., 2007;Kussakin, 1973;Strugnell et al., 2008). Contrastingly, for some deep-sea taxa, such as protobranch bivalves, gastropods, and some crustaceans, depthrelated diversification have been observed, indicating possible depth-related barriers to dispersal (Etter et al., 2011Etter and Bower, 2015;Havermans et al., 2013;Wilson, 1983). The depth-related population differentiation observed in these studies, however, covered larger bathymetric ranges than the ones studied here and were mostly situated at the lower bathyal and abyss (Etter et al., 2011;Etter and Bower, 2015). Just a few studies have assessed shallow-deep connectivity using a combined morphological and molecular approach (Bik et al., 2010;Riehl and Kaiser, 2012;Van Campenhout et al., 2014;Van Gaever et al., 2009). In contrast to the molluscs and crustaceans mentioned above, selected nematodes and isopods show high degrees of genetic similarity across depths, suggesting taxon-specific barriers (Bik et al., 2010;Riehl and Kaiser, 2012).
Except for Halalaimus, all deep-sea nematode genera sequences obtained in this study were sequenced for the first time. Although the relatively conserved 18S rDNA used here may not be the most suitable marker to assess dispersal, evolutionary rates of this gene are unknown for the nematode genera studied. Nevertheless, the presence of identical sequences between individuals from shallow and deep habitats (Fig. 8) provides hints towards dispersal between depths at relatively recent evolutionary timescales. Our results for Halalaimus are in accordance with Bik et al. (2010), revealing multiple historic interchanges between habitats of different depth for multiple species. Likewise, no clear geographical structuring was observed in our phylogenetic tree, although this result could be biased due to limited taxon and geographic sampling. Whether nematode dispersal occurs passively through hydrodynamics or is active employing chemical cues and active swimming, exchange among marine nematode assemblages can be maintained both over large (> 500 km) and small (50-100 km) geographical distances. This explains the success of these benthic organisms as colonizers (Boeckner et al., 2009;Derycke et al., 2013;Gallucci et al., 2008).
Even though our results indicate a link between shallow and deep habitats, other studies have suggested endemism in deep-sea habitats for nematodes (De Mesel et al., 2006;Van Campenhout et al., 2014). For example, Halomonhystera disjuncta was previously believed to occur in both shallow and deep habitats (Van Gaever et al., 2009), but a recent study based on 18S, COI, and ITS sequences showed that this species in fact constitutes two different lineages occupying deep and shallow environments, respectively (Van Campenhout et al., 2014). The bathymetric and geographic range was, however, broader than here.
Our phylogenetic results moreover highlight the scarcity of publicly available DNA sequence data for deep-sea nematodes. For example, we present here the first sequence of the genus Microlaimus (no records in the GenBank, searched on 3 December 2015). Other genera are poorly represented in public sequence depositories (e.g. Gammanema: two 18S sequences; Leptolaimus: three 18S sequences; and Richtersia: two 28S sequences). The use of more specific and variable markers, such as the mitochondrial COI or the rDNA internal transcribed spacer, was not possible within this study due to low success rate of DNA amplification. Low success rates in PCR amplification are a known issue in deep-sea nematodes, but the causes are not well understood (Bik et al., 2010). Degradation of DNA may have occurred during sample processing and could be caused by increases in temperature.
Genetic structuring of shallow-water nematode populations was shown by Derycke et al. (2013) based on more variable markers (COI, ITS). They showed that despite being capable of long distance dispersal, nematodes may also show clear genetic differentiation at small spatial scales. In this study, we demonstrated the high dispersal capabilities and connectivity for nematodes, but those were not high enough to counteract community differentiation observed in the genus composition. Moreover, it is likely that only a small number of species show relatively high dispersal (gene flow), while other species may have limited dispersal abilities.

Conclusions
Sediment heterogeneity seemed to be the most important factor responsible for the greater variation in nematode community structure within and between stations from the same transect, as well as between transects. "Deep" stations were uniformly dominated by silt clay and exhibited low nematode turnover and less variable nematode communities typically found in deep-sea soft sediments, dominated by the genera Acantholaimus and Halalaimus. High sediment variability at the "shallow" stations was associated with a more diverse and variable nematode community between stations, where high turnover was characterized by the dominance of different genera at each station. This high sediment heterogeneity at the shelf break suggests that strong near-bottom current pulses might be an indirect factor promoting diversity at the scale of transect through the creation of patches and redistribution of resources (Fig. 9). Moreover, the high beta diversity observed between transects surpassed the turnover within the shelf break and upper slope transects. These differences were attributed to the even higher differences in sediment variability observed between these two transects. This effect of sediment variability on the nematode community suggests that environmental factors were the main responsible for differences observed between stations and transects, rather than geographical distance or depth (Fig. 9). Nevertheless, despite the larger differences observed between transects, still a considerable percentage of shared genera were observed between the shelf break and the upper slope. Finally, phylogenetic clusters suggested that depth may not be a factor isolating populations of the nematode genus Halalaimus.

Data availability
The DNA sequence data is publicly accessible in GenBank (Benson et al., 2008).
The Supplement related to this article is available online at doi:10.5194/bg-14-651-2017-supplement.