Sediment Features and Human Activities Structure the Surface Microbial Communities of the Venice Lagoon

Sediment microbial communities play essential roles in marine ecosystem functioning. Their study is crucial to understand how environmental conditions affect microbial diversity and biogeochemical cycles. The Venice Lagoon, one of the largest Adriatic transitional systems, is subjected to different natural and anthropogenic stressors. In this study, surface sediments microbial communities were finely characterized using 16S rRNA gene amplicon sequencing and metagenomics. Samples were collected seasonally in 2019 and 2020 at different sites within Lagoon sub-basins. Our results indicated a stable spatial distribution of the sediment microbiome: salinity, grain size, and total organic carbon were found to be important drivers in shaping prokaryotic communities, while temperature had a minor role. We detected high microbial diversity at all stations, mainly due to low abundant taxa: bacteria represented the vast majority of the reads (∼96.1%), with Woeseia the most abundant genus (∼4.4%). The metagenomic analysis highlighted significant differences among sites in terms of biogeochemical processes (e.g., C, N, Fe, and S metabolism), and cell-cell interaction strategies (e.g., mobilome, regulations and cell signaling). Chioggia, a sandy site subjected to marine influence, presented the highest abundance of ammonia-oxidizing archaeon Candidatus Nitrosopumilus, in accordance with the highest amount of ammonia monooxygenase subunit genes. At the same site, sulfate-reducing bacteria (Desulfobacteria and Desulfobacterales) and sulfur-related genes were found in lower abundance. Marghera and Tresse, the most polluted sites, showed higher abundance of sewage-related bacteria and antibiotic and toxic compound resistance genes. Furthermore, these sites showed higher amount of genes related to cell-cell interaction, such as pathogenicity islands, transposable-elements, and biofilm formation. Our findings highlighted that sediment features and human-related activities have profound and long-term impacts on the surface sediment microbial communities of the Venice Lagoon.

Sediment microbial communities play essential roles in marine ecosystem functioning. Their study is crucial to understand how environmental conditions affect microbial diversity and biogeochemical cycles. The Venice Lagoon, one of the largest Adriatic transitional systems, is subjected to different natural and anthropogenic stressors. In this study, surface sediments microbial communities were finely characterized using 16S rRNA gene amplicon sequencing and metagenomics. Samples were collected seasonally in 2019 and 2020 at different sites within Lagoon sub-basins. Our results indicated a stable spatial distribution of the sediment microbiome: salinity, grain size, and total organic carbon were found to be important drivers in shaping prokaryotic communities, while temperature had a minor role. We detected high microbial diversity at all stations, mainly due to low abundant taxa: bacteria represented the vast majority of the reads (∼96.1%), with Woeseia the most abundant genus (∼4.4%). The metagenomic analysis highlighted significant differences among sites in terms of biogeochemical processes (e.g., C, N, Fe, and S metabolism), and cell-cell interaction strategies (e.g., mobilome, regulations and cell signaling). Chioggia, a sandy site subjected to marine influence, presented the highest abundance of ammoniaoxidizing archaeon Candidatus Nitrosopumilus, in accordance with the highest amount of ammonia monooxygenase subunit genes. At the same site, sulfate-reducing bacteria (Desulfobacteria and Desulfobacterales) and sulfur-related genes were found in lower abundance. Marghera and Tresse, the most polluted sites, showed higher abundance of sewage-related bacteria and antibiotic and toxic compound resistance genes. Furthermore, these sites showed higher amount of genes related to cellcell interaction, such as pathogenicity islands, transposable-elements, and biofilm formation. Our findings highlighted that sediment features and human-related activities have profound and long-term impacts on the surface sediment microbial communities of the Venice Lagoon.

INTRODUCTION
Globally, marine microbes living in the top 10-50 cm account for 1.7-5.0 × 10 28 cells (Whitman et al., 1998;Danovaro et al., 2015;Flemming and Wuertz, 2019). Within transitional shallowwater benthic ecosystems, the surface sediment is a highly spatially and temporally dynamic habitat. It is subjected to frequent resuspension and deposition events, mixing, transport, and it is shaped by significant biotic activities of macro-to microorganisms (Daumas, 1990;Mayer, 1993). The surface sediment is highly heterogeneous and characterized by high water content, steep physical and chemical gradients that structure the ecology, the interactions, and the metabolism of microorganisms (e.g., light, oxygen, and organic matter) (Zinger et al., 2011). The microbes within this layer play an essential role in the decomposition of organic matter, nutrient cycling, benthic food webs, thus governing the ecosystem functioning (Schallenberg and Kalff, 1993;Azam and Malfatti, 2007). In respect to the microbial diversity and functions of the pelagic realm, fewer information is available on benthic microbial communities (Probandt et al., 2017). Surface sediments usually present higher microbial abundance and diversity in terms of taxa, potential functions, and metabolisms (Zinger et al., 2011) in comparison to the deeper ones. This microbial environment is influenced by the water column biogeochemistry, and on the other side, from the deeper sediment biogeochemistry, via porewater fluxes (e.g., anoxic conditions, reduced state of organic matter, and contaminants). The study of the microbial distribution patterns within "sediment skin" is fundamental to understand and predict the responses of the marine ecosystem to environmental changes (Zinger et al., 2011). In contrast to the water environment (with the exception of marine snow particles), microbial abundances reach 10 9 cells for cm 3 of sediment and within these universes there is harsh microscale competition for resources (e.g., energy and nutrients) and space (Petro et al., 2017).
Microbes exert their ecosystem role at the microscale (Azam and Malfatti, 2007;Stocker et al., 2008) displaying a great variety of adaptive strategies from dormancy to exchanging electrons and specialized behaviors. Specifically, motility, twocomponent systems, antiviral systems, the propensity to exchange genetic material, secretion system, and production of redox active antibiotics may promote overall their persistence and success in such crowded environment.
The Venice Lagoon (northern Adriatic Sea), with a surface of 550 km 2 , is one of the largest transitional systems in the Mediterranean Sea (Madricardo et al., 2019). This lagoon includes different environments, such as salt marshes, channels, fish farms, natural and artificial islands. It is subjected to different and extensive natural and anthropogenic stressors (Solidoro et al., 2010) such as sea-level rising and subsidence, elevated touristic pressure on Venice city (Seraphin et al., 2018), and the presence of industrial-chemical areas (Bellucci et al., 2002). The essence to be a lagoon relies on the tight coupling between the water and sediments; the sediment environment provides information on the dynamic history of these transitional areas. In the Venice Lagoon, several contaminants have been detected, coming mainly from industrial processes, incinerators, drainages, and urban wastes (Solidoro et al., 2010). These include heavy metals, toxins, polychlorinated biphenyls (PCBs), and others (Han et al., 2007;Borin et al., 2009;Solidoro et al., 2010;Gieskes et al., 2015). These pollutants were discharged in the water column and then they accumulated in the sediments, by the interaction with the downward fluxes of particulate matter (DePinto et al., 2010). The quality of the Venice Lagoon sediments has been evaluated via different risk assessment screenings: several contaminants exceed sediment quality guidelines thresholds, with the highest screening risk level in the north and central part of the lagoon, near the two hotspots of contamination Porto Marghera and Venice city canals (Apitz et al., 2007).
Even if the Venice Lagoon is one of the most studied coastal ecosystems (Solidoro et al., 2004), more effort needs to be done by employing molecular approaches focusing on the microbial role within this transitional area. So far, a handful of high-profile studies have been carried out, but none of them have employed the metagenomic approach. In 2007, Danovaro and Pusceddu used automated ribosomal intergenic spacer analysis (ARISA) to study microbial diversity in some Mediterranean coastal lagoons, including Venice. They found no relationships between the environmental features of the lagoon and the bacterial diversity, but, on the other hand, the latter was significantly and positively correlated with the ecosystem functioning and efficiency (Danovaro and Pusceddu, 2007). Borin et al. (2009) used 16S rRNA gene clone libraries to assess microbial community structure and diversity in the different subbasins of the Lagoon. They highlighted that the sediments, mainly anoxic, were colonized by microbial communities that presented a species richness correlated with total elemental sulfur and acidvolatile sulfide (Borin et al., 2009). More recently, Quero et al. (2017) applied DNA metabarcoding to Venice Lagoon planktonic and benthic bacterial assemblages to explore diversity patterns and the role of environmental gradients. They found evidence of different temporal dynamics in the pelagic and benthic domains (Quero et al., 2017).
This study aims to provide new insights into the spatial and temporal distribution and dynamics of microbes living in the surface sediment of the Venice Lagoon in terms of diversity, functions, and metabolisms. The metagenomic approach was employed for the first time in the Lagoon, enhancing our knowledge of the behavior of the surface sediment microbial communities. Our data contribute to the understanding of the importance of the surface sediment within the microbial biogeochemical cycles of carbon, nutrients, metals, pollutants, antibiotics, and xenobiotics from the large to the micro-scale. Moreover, we have used the sediment microenvironment to test hypotheses on the microscale functional adaptations of microbes, shedding light on the potential behavior in such a complex and heterogeneous matrix.

Sampling
The sampling was carried out seasonally in 2019 (Supplementary Table 1) at five sites of the Venice lagoon, Italy (Figure 1): Frontiers in Marine Science | www.frontiersin.org Chioggia (C), Marghera (M), Palude della Rosa, (P), Sacca Sessola (S), and Tresse (T). Furthermore, C, M, and P were sampled with the same time frequency in 2020. The sites were chosen to cover the four sub-basins, defined according to risk analyses following the international sediment quality guideline of the Lagoon (Apitz et al., 2007;Supplementary Table 2). They are located in the southern (C), in the central (S), in the centralnorth (M and T), and in the northern (P) sub-basins. The Chioggia site is near the city of Chioggia, very close to the southern lagoon inlet and more subjected to marine influence. The Sacca Sessola site, situated in the center of the Lagoon, south of Venice city, is potentially affected by urban contamination (Picone et al., 2016). The Marghera site is inside the active industrial area of Porto Marghera and is characterized by high anthropogenic impact. The Tresse site is next to an artificial island created as a dumping site for urban and industrial waste. Finally, Palude della Rosa site is less impacted by human activities (Picone et al., 2016).
In Chioggia, Marghera, and Palude della Rosa samples were collected with a Van Veen grab (0.045 m 2 ), while in Sacca Sessola and Tresse with manually operated corers (10 cm diameter). The first centimeter of sediment was collected with a sterile spatula. The sampling was performed in three independent replicates. Water temperature and salinity were recorded with the multiparametric probes HI98128 (HANNA Instruments, Italy) and AP-2000 (Aquasearch, Italy) (Supplementary Table 3).

Sediment Characterization and Bacterial Abundance Measurement
The water content was estimated by weighing fresh and dried (105 • C for 24 h) sediment aliquots (2 cm 3 ). The grain size composition was determined with a laser Particle Size Analyzer (BECKMAN COULTER LS 13 320) after treating 1 g of wet sediment with H 2 O 2 . The sediment was described using Shepard's classification (Shepard, 1954) and the data were expressed as percentages of sand, silt, and clay. The limit between silt and clay was fixed to 4 µm, following the Udden-Wentworth classification (Wentworth, 1922).
For Total Organic Carbon (TOC) estimation, freeze-dried sediment was grounded in a ceramic mortar and sieved on a 250 µm iron steel sieve (Endecotts LTD., United Kingdom). Triplicate subsamples of about 8-12 mg were weighed on a micro ultrabalance with an accuracy of 0.1 µg. Before TOC determination, subsamples were treated directly into capsules with increasing concentrations of HCl (0.1 N and 1 N) to remove the carbonate fraction (Nieuwenhuize et al., 1994). Carbon content was determined using a CHNO-S elemental analyzer ECS 4010, Costech, Italy) according to Pella and Colombo (1973). Standard acetanilide (Costech, purity ≥ 99.5%) was used to calibrate the instrument and empty capsules were also analyzed to correct for blank. Measurement quality control was performed using internal standards and it was also verified for carbon against the certified marine sediment reference material PACS-2 (National Research Council Canada). Significant differences in TOC were calculated in the R environment (v. 4.0.3, R Core Team, 2019) by Kruskal-Wallis and Wilcoxon non-paired tests.
False discovery rate (FDR) was used for p-value correction. Results with q-value < 0.05 were considered significant.
The prokaryotic abundance was estimated by flow cytometry following Deng et al. (2019) with modifications, using a FACSCanto II (Becton Dickinson) instrument equipped with an air-cooled laser at 488 nm and standard filter setup. A sediment slurry was subjected to mechanical shaking, ultrasonication, and hydrofluoridric acid (1%, final concentration) treatment. Slurry mixed with the Stop solution was filtered through a 10 µm polycarbonate filter, diluted, and stained with SYBRGreen I (2×, final concentration). Negative controls were prepared by autoclaving sediments and running unstained sediments. The acquisition threshold was set to green fluorescence and stained prokaryotic cells were identified in the side scatter versus green fluorescence plot. The flow rate was calibrated daily, by running distilled water and weighing it before and after the run (at least five replicates). Data were acquired and processed with the FACSDiva software (Becton Dickinson). Prokaryotic cells per gram of dry sediment were calculated using the acquired cell counts, the flow rates and the respective sample water content, and statistical analyses were run as described for TOC data.
A Principal component analysis (PCA) was performed to explore and visualize similarities among the samples based on environmental variables (temperature, salinity, grain size, TOC and prokaryotic abundance) with the R package FactoMineR (Le et al., 2008), factoextra (Kassambara, 2015), and ggplot2 (Wickham et al., 2016).

DNA Extraction and Sequencing
DNA was extracted with DNeasy PowerSoil Pro Kit (Qiagen) following the manufacturer's instructions. Quality and quantity of the extracted DNA was assessed with Nanodrop spectrometer (Thermo Fisher Scientific) and with Qubit Fluorimeter (Thermo Fisher Scientific).
For the amplicon sequencing analysis, sediments from both 2019 (at each site, one sample per season for a total of 20, in three replicates) and 2020 (at C, M, and P sites, one sample per season for a total of 12, in three replicates) were considered. The V4-V5 region of 16S rRNA gene was amplified using 515-Y (5 -GTGYCAGCMGCCGCGGTAA-3 ) and 926R (5 -CGYCAATTYMTTTRAGTTT-3 ) primers (Parada et al., 2016). Libraries were prepared following the 16S Metagenomic Sequencing Library Preparation protocol and run on an Illumina MiSeq System for a read length of 2 × 250 bp at the genetic and epigenetic ARGO Open Lab Platform, Area Science Park, Trieste, Italy.
For the metagenomic analysis, sediments from 2019 (at each site, one sample per season for a total of 20, in three replicates) were considered. Libraries were prepared following the Illumina Nextera DNA Flex Library Prep protocol and run on an Illumina NovaSeq 6000 System for a read length of 2 × 300 bp at the genetic and epigenetic ARGO Open Lab Platform, Area Science Park, Trieste, Italy.

16S Bioinformatic Pipeline
Bioinformatic analyses were performed with QIIME2 2020.6 (Bolyen et al., 2019). Raw sequences were quality filtered and FIGURE 1 | Sampling sites in the Venice Lagoon. The maps were created with the software Ocean Data View (Schlitzer, 2018). denoised with DADA2 (Callahan et al., 2016). Alpha-diversity metrics were estimated after samples were rarefied. To assess the reliability and robustness of the replicates, a correlation analysis (Spearman's ρ) on the amplicon sequence variants (ASVs) table was performed with the R package ggplot2 (Wickham et al., 2016). Taxonomy was assigned to ASVs using the sklearn naïve Bayes taxonomy classifier (Bokulich et al., 2018) against the Silva 138 99% reference database with 7-level taxonomy (Quast et al., 2013). Reads belonging to Eukarya, mitochondria, chloroplast, and with frequency < 2 (singletons) were removed.
Statistical analyses were conducted in the R environment (v. 4.0.3, R Core Team, 2019). Rarefaction curves on ASVs table were plotted with the R package vegan (Oksanen et al., 2019). Significant differences in alpha diversity metrics were calculated by Kruskal-Wallis H test. FDR was used for p-value correction. Results with q-value < 0.05 were considered significant.
To visualize similarity patterns of prokaryotic communities, a Principal Coordinates Analysis (PcoA) using Bray-Curtis dissimilarity matrices was constructed using normalized read abundances in the R package vegan (Oksanen et al., 2019). A permutational multivariate analysis of variance (PERMANOVA) with 4999 permutations was computed on normalized ASVs table to investigate the effect of site and season using the function adonis in the R package vegan (Oksanen et al., 2019).
To further explore the environmental drivers of the community, we used a distance-based redundancy analysis (dbRDA) using the capscale function in the R package vegan (Oksanen et al., 2019). The environmental variables (temperature, salinity, grain size, and TOC) were normalized using a z-score transformation and variance partitioning analysis was used to quantify the contribution made by each variable using the varpart function in the R package vegan (Oksanen et al., 2019). The significance of the variables was tested with ANOVA and 4999 permutations after building a reduced model. As the grain size components (sand, silt, and clay) present multicollinearity, only one component (sand) was included in the analysis.
An indicator species analysis was performed at site level with the R package indicspecies (De Cáceres and Legendre, 2009). The results were visualized as networks by mean of the R package igraph (Csard and Nepusz, 2006). Moreover, an analysis of composition of microbiomes (ANCOM; Mandal et al., 2015) was used to test for differences at the genus level. Microbial diversity visualization at different taxonomic level was done using the average normalized abundance of the three replicates of each sample, with R package phyloseq (McMurdie and Holmes, 2012) and ggplot2 (Wickham et al., 2016). Microbial signatures of fecal-and sewage-associated bacteria (McLellan et al., 2010;Newton et al., 2013;Paliaga et al., 2017) were also evaluated. Statistical Analyses of Metagenomic Profiles 2.1.3 (STAMP; Parks and Beiko, 2010) was used to assess differences in the relative proportion of sequences belonging to the different taxa considering the whole dataset. Significant differences were calculated by ANOVA and Tukey-Kramer post hoc test. FDR was used for p-value correction. Results with q-value < 0.001 were considered significant.

Metagenomic Analysis
The number and quality of reads were checked with Fast QC (Andrews, 2010). Reads were assembled into contigs using MEGAHIT-1.2.9 (Li et al., 2015) combining the three replicates of each sample after testing the significance of their correlations. Gene prediction was performed with Prodigal 2.6.3 (Hyatt et al., 2010) on contigs longer than 1,000 bp ( Van der Walt et al., 2017).
The contigs > 1,000 bp were submitted to Metagenome Analysis Rapid Annotation using Subsystem Technology (MG-RAST; Keegan et al., 2016) for a similarity search using the SEED Subsystem (Overbeek et al., 2005) and the Kyoto Encyclopedia of Genes and Genomes (KEGG; Kanehisa and Goto, 2000) databases, using a cut-off E-value of 1e −5 , minimum identity of 60% and a minimum alignment length of 15 bp. SEED and KEGG annotations were normalized by the number of reads assigned to the prokaryotic single copy gene recA (Acinas et al., 2021) prior to calculating the relative abundances.
Key oceanic marker genes (Sunagawa et al., 2015) and key sediment metabolic genes (Dombrowski et al., 2018) were evaluated to provide an overview of the main biochemical cycling metabolisms present. We have created a novel ad hoc specific microscale marker function list to explore the microbial potential at the small-scale in the sediment environment.
Annotation results were visualized using Statistical Analyses of Metagenomic Profiles 2.1.3 (STAMP; Parks and Beiko, 2010) to assess differences in their relative proportion considering the whole dataset. Significant differences among samples were calculated by ANOVA and Tukey-Kramer post hoc test. FDR was used for p-value correction. Results with q-value < 0.001 were considered significant.

Sediment Characterization
The five stations were characterized by different percentages of the grain size classes of sand, silt, and clay (Supplementary Table 4), in agreement with previous studies (Borin et al., 2009;Picone et al., 2016;Zonta et al., 2018). Based on Shepard's classification, Chioggia was categorized as sand, Marghera as clayey silt, Palude della Rosa as silty-loam, Sacca Sessola as loam, and Tresse as silty sand.

Prokaryotic Diversity and Community Composition
A total of 31,189,504 raw sequences were generated. After the denoising procedure, 20,191,333 were retained with an average of 219,471 ± 92,086 sequences per sample. The average number of ASVs was 5,230 ± 1,536. The rarefaction curves indicated that the sequencing effort was enough to assess the sample biodiversity (Supplementary Figure 2). The replicates significantly correlated (p < 0.05) with an average ρ of 0.86 ± 0.09.
Shannon's diversity index showed an average of 10.7 ± 0.4 and Pielou's evenness showed an average of 0.89 ± 0.02 (Supplementary Figure 3). Overall, Autumn had a significantly higher index than the other seasons. No significant differences in biodiversity indices were detected among sites and siteseason combinations.
The PERMANOVA highlighted that site, season, and their interaction had a significant effect (p < 0.001) on the prokaryotic community structure (Supplementary Table 6). The highest amount of variance was explained by the site (37%), followed by season-site interaction (11%) and season (6%).
The PCoA (Supplementary Figure 4), in agreement with the PERMANOVA, grouped the samples by site.
The variance partitioning analysis showed that environmental parameters explained ∼20% of the variance (Supplementary Table 7). The dbRDA (Figure 3) revealed that grain size, TOC, and salinity had a significant role (p < 0.001) in shaping the communities, while temperature was not a significant variable (p > 0.05).
Regarding the taxonomic composition, at the domain level most reads were assigned to Bacteria (average on the overall dataset 96.1 ± 2.6%), and a smaller percentage (3.8 ± 2.6%) to Archaea.
Among the genera related to microscale processes, we detected the cable bacteria Candidatus Electrothrix at all sites (0.1 ± 0.2%); these bacteria have long, multicellular filaments that can conduct electric currents over centimeterscale distances (Trojan et al., 2016). Likewise, lignin-degrader genera (Bacillus, Novosphingobium, Pseudomonas, Rhodococcus, and Streptomyces; Lee et al., 2019) were found at all sites (0.006 ± 0.008%). The Indicator Species analysis showed that each site hosted a set of indicator bacterial genera ( Figure 5, Supplementary Table 8), with Chioggia showing the highest number (n = 21). On the contrary, excluding Chioggia, the highest number of genera (n = 27) were shared among Marghera, Palude della Rosa, Sacca Sessola and Tresse. Chioggia was characterized by the presence of aerobic (or strictly aerobic), halophilic bacteria such as Aquimarina, Aureisphaera, and Pseudoalteromonas and by the absence of genera of sulfate-reducing Desulfobacteraceae that were shared among the other sites, like Desulfoconvexum, Desulfonema, and Desulfobacter. ANCOM analysis (Supplementary Figure 6) confirmed the significant higher abundance ammonia-oxidizing archaeon Candidatus Nitrosopomilus in Chioggia and a lower presence of the Desulfobacteraceae member Desulfobacter. Moreover, it highlighted a higher presence of the nitrite-oxidizing Nitrospira and of the aerobic Ahrensia and Granulosicoccus. Palude della Rosa showed significantly higher presence of the halophile Halochromatium and the glycogen-accumulating Candidatus Competibacter.
Microbial signatures of selected anthropogenic-related taxa (Figure 6) presented overall 0.1-0.3% relative abundance for fecal-associated bacteria and 0.01-0.04% for sewage-associated bacteria. In Sacca Sessola and Chioggia, fecal-associated microbes showed the highest presence, mainly due to Clostridiaceae, while Marghera showed the highest proportion of Bacteroidaceae and Ruminococcaceae. For the sewage-associated bacteria, Marghera showed on average the highest proportion of reads due to Acinetobacter and Trichococcus while Palude della Rosa showed the lowest proportion.

Metagenomic Analyses and Functional Characterization
A total ∼120 million of raw reads (30 Gbp) were generated, specifically 14.8 ± 3.7 millions of paired-end reads with quality score of 35 ± 2 for each metagenome. On average, contigs derived from the three co-assembled replicates for each sample were 5,208,341 ± 376,210 contigs > 1,000 bp, and 652,755 ± 174,565 coding sequences.
Key oceanic marker genes (Sunagawa et al., 2015) and key sediment metabolic genes (Dombrowski et al., 2018) have been selected and clustered as follows: Carbon fixation, Energy and metabolism, and Carbon-nutrient coupling.
Regarding anthropogenic-related functions, the "Virulence, Disease and Defense" category was significantly higher in Marghera and Tresse in respect to the other sites (3.36 ± 0.02% vs. 3.04 ± 0.03%), mainly due to a significantly higher presence of genes in the "Resistance to antibiotics and toxic compounds" group (3.08 ± 0.04% vs. 2.72 ± 0.03%).
The key specific microscale marker functions, identified with the SEED database, have been clustered as follows: DNA metabolism (CRISPs and DNA uptake and competence), Mobilome (Gene Transfer Agents, integrons, pathogenicity islands, phages and prophages, plasmids, and transposable elements), Regulation and Cell Signaling (biofilm formation, phenazine biosynthesis, osmotic stress, sporulation, siderophores, toxic-antitoxic systems, and two-component systems) and Nanomachines (protein secretion system Type IV and Type VI) (Figure 9).

Microbial Community Diversity and Structure Are Shaped by Site Features
Microbial communities in marine sediment are extremely diverse due to the combination of different factors such as habitat temporal stability, high niche diversity, and resource partitioning due to the complex physical-chemical gradients (Zinger et al., 2011;Acosta-Gonzalez and Marques, 2016). As such, they represent a major reservoir of genetic variability in the marine environment (Polymenakou et al., 2005).
Our analysis detected a pronounced diversity in the communities, in agreement with 16S studies in the Venice Lagoon (Borin et al., 2009;Quero et al., 2017) and in other transitional and coastal benthic systems (Bianchelli et al., 2020;Miksch et al., 2021). The sampling site specific features (grain size, salinity and TOC content, Figure 2) were important drivers in shaping the prokaryotic assemblages, while seasonality had a significant but minor role (Figure 3). Depending on the area, seasonality in the sediment microbial communities is not as regularly and clearly found as in the water column (Tšertova et al., 2011;Gobet et al., 2012;Miksch et al., 2021), even in the Venice Lagoon (Quero et al., 2017). Quero et al. (2017), in a study including pelagic and benthic communities in the central part of the Venice Lagoon, found seasonality in both compartments but less pronounced in surface sediments than in the overlaying water, where this feature is well documented (Celussi et al., 2009).
In marine sediments, it is now known that physical-chemical characteristics including salinity, chlorophyll a, TOC content, grain size, and benthic fauna can influence the microbial communities (Sapp et al., 2010;Pala et al., 2018;Hoshino et al., 2020). Generally, the composition of the sediment bacterial communities can present high dissimilarities between sites, reflecting the heterogeneity of these environments, mostly due to limited physical mixing and the geographical distance (Acosta-Gonzalez and Marques, 2016).
We found that salinity, grain size, and TOC significantly contributed to the structure of the microbial communities in the different sites of the Venice Lagoon (Figure 3).
Salinity is recognized as a major environmental determinant of microbial composition in benthic communities (Lozupone and Knight, 2007;Bolhuis et al., 2013) in coastal (Bolhuis and Stal, 2011;Severin et al., 2012;Bolhuis et al., 2013) and lagoon systems (Tsuboi et al., 2013;Pavloudi et al., 2016;Behera et al., 2017). In our study, the more marine-influenced site (Chioggia, located at the southern inlet of the Venice Lagoon) was confirmed as the most diverse one, forming a separated group in the multivariate analyses (see Figure 3).
It is known that grain size can affect microbial abundance and activity by providing different surface properties (e.g., fine sediments provide greater areas that can be colonized) and influencing permeability, water chemistry, hydraulic conductivity, and pore space between sediment particles (Dale, 1974;Zhang et al., 1998). A positive correlation between bacterial abundance, sediment grain size, and organic carbon content is well established (Dale, 1974) and reported in different studies (Santmire and Leff, 2007;Legg et al., 2012;Tsuboi et al., 2013;Fazi et al., 2020). In another Italian Adriatic lagoon, it was found that lowest salinity, higher percentage of fine grain sediment and organic matter corresponded to the highest microbial abundances (Cibic et al., 2019). This was confirmed also in our study, in which the sites enriched in silt, Marghera and Palude della Rosa, presented the highest prokaryotic abundance and TOC concentration.
Sediment grain type influences indirectly the microbial assemblage at the microscale, due to the variable water fluxes (Ahmerkamp et al., 2020) that changes the redox conditions thus driving more diverse and able-to-adapt communities in sandy sites compared to the clay and silty ones.
Being and important source of energy for heterotrophs the organic carbon (TOC) pool is an important structuring agent for benthic microbial communities (Oni et al., 2015). As such, in our study the content of TOC was important in structuring the differences among the communities (e.g., Chioggia vs. other sites), even though we did not analyze its quality (biopolymers' concentration, labile vs. recalcitrant matrices, etc.), which is also known to affect the composition of bacterial communities (e.g.,   Fazi et al., 2020). Likewise, ecological factors that were not investigated in this study can contribute to the differentiation of the sites: interspecific bacteria competition, grazing pressure, viral lysis, vegetation, contaminants, and also stochastic events (Böer et al., 2009;Luna et al., 2013;Acosta-Gonzalez and Marques, 2016;Bianchelli et al., 2020).

The Rare Biosphere Supports Microbial Diversity
The microbial communities were characterized by few dominant taxonomic groups, the core microbiome, and a high number of low abundant taxa, the so-called rare biosphere (Sogin et al., 2006), that accounted for most of the phylogenetic diversity.
For instance, among the 84 phyla detected in the whole dataset, 49 were considered rare, as represented by an average relative abundance < 0.01% (Galand et al., 2009;Gobet et al., 2012), as well as for 1181 out of the 1587 detected genera. The presence of such rare biosphere is of ecological importance as these taxa, often overlooked, can have an over-proportional role in biogeochemical cycles and can represent a hidden driver of ecosystem functions (Jousset et al., 2017).
Even if we detected significant differences among sites in high abundant taxa, from the phylum to the genus level, the numerous low-abundant taxonomic units were those that more contributed to their site-specific fingerprint (see Figures 3, 5). These taxa were mostly uncharacterized (e.g., belong to uncultured bacteria typically detected only with molecular techniques and are still not defined in terms of physiology and metabolism), and more studies are needed to define their ecological roles.
In accordance with Jousset et al. (2017) we suggest that the rare biosphere in the Venice Lagoon can be considered as an ensemble of keystone species that provide stability and seeds for ecosystem functions, specifically in terms of water purification, nutrients, and global climate regulations (da Mosto et al., 2020).
Regarding the most abundant taxa, the globally distributed (Gamma)Proteobacteria were confirmed to be prevalent at all sites (see Figure 4), as generally detected in sediment in both the costal and the deep seafloor (Fry et al., 2008;Zinger et al., 2011;Acosta-Gonzalez and Marques, 2016), including the Venice Lagoon one (Borin et al., 2009;Quero et al., 2017). Among these, the ubiquitous genus Woeseia was the most abundant (see Figure 4). Recent molecular characterization suggested that it is likely to grow on proteinaceous matter, potentially derived from detrital cell membranes, cell walls, and other organic remnants Hoffmann et al., 2020). Several significant differences in taxonomic composition were highlighted by ANOVA/Tukey-Kramer and ANCOM analyses in sulfur-related bacteria (Desulfobacteria, Desulfobacterales) involved in sulfate reduction. In our study, these taxa were less represented in Chioggia where, differently form the inner area of the Lagoon where the sediment is mostly anoxic, hypoxic conditions are commonly found (Borin et al., 2009). Chioggia sandy sediment is in fact more exposed to the effect of winds, tidal currents, and ship movements, that favor partial oxygenation (Borin et al., 2009). Permeable sandy sediment pore space is constantly flushed by the overlying water and could trap detritus and living cells from the water column (Boudreau et al., 2001;Gobet et al., 2012). In the other sites, the higher amount of organic matter contributes to the oxygen depletion and to the activity of sulfate-reducing bacteria in sulfide production (Zaggia et al., 2007). This is further supported by the metagenomics analysis, which also highlighted a significantly lower number of genes related to sulfur metabolism in Chioggia (see section "Wide Range of Carbon Fixation and Carbon Flow Strategies").

Wide Range of Carbon Fixation and Carbon Flow Strategies
Sediment bacteria play a key ecological and biogeochemical role in marine ecosystems. They present high abundance and genetic variability, critical role in the transformation and speciation of major bioactive elements (e.g., carbon, nitrogen, phosphorus, oxygen, and sulfur), and in the degradation of organic pollutants (Polymenakou et al., 2005;Wu et al., 2008). The metagenomics approach adopted in this study allowed the identification of the genes involved in biogeochemical cycles providing insights into the potential metabolic functioning of the microbial communities (Doney et al., 2004;Dombrowski et al., 2018;Acinas et al., 2021). Following the "carbon currency" in the microbial sediment compartment we were able to map the carbon fixation processes, the energy and metabolism processes, the carbon-nutrient coupling, the carbon reworking processes, which includes the degradative processes, and the behavior and anthropogenicrelated metabolisms (see Figure 7).
In our sites, we found that oxygenic photosynthesis (by Prokaryotes) was the main carbon fixation pathway. Anoxygenic photosynthesis was also detected together with five additional carbon fixation pathways which differ in reducing compounds, energy source, and oxygen sensitivity of enzymes alternative pathways (Wood-Ljungdahl, Hydroxypropionate bicycle, hydroxypropionate/hydroxybutyrate, C1-metabolism). These findings highlighted the metabolic flexibility of the community in fixing CO 2 from different sources thus suggesting high diversity in C fixation modes in the sediment. Once the carbon is turned into biomass and organic carbon, it can be metabolized via aerobic and anaerobic respiration or fermentation. The surface sediment is affected by the light diurnal cycle and by the presence of oxygen (as a by-product of photosynthesis and from the oxygenated seawater column) (Petersen et al., 1994;Nealson, 1997). Aerobic respiration was the main pathway to respire the organic matter and secondly, anaerobic respiration (N, Mg, Fe, S; i.e., changing the terminal electron acceptor according to their availably, Nealson, 1997) was detected. Furthermore, often understudied metabolisms, such as methanogenesis, butane, proton degradation (hydrogenase, with a great diversity of enzymes), and fermentation were also important in enriching the metabolic repertoire of the microbial communities across sites (Falkowski et al., 2008). Within the concept of shared metabolisms, cable bacteria (family Desulfobulbaceae, Kjeldsen et al., 2019) were present in our samples. These microbes are formed by centimeter long filaments that span from the oxic water-sediment interface to the sulphidic deeper sediment layer. Within their filament, they couple sulfide oxidation with oxygen or nitrate reduction via long-distance electron transport.
Sulfur-related genes were found in high abundance in our sites with the exception of the sandy site of Chioggia. This is possibly due to the sediment texture that is more permeable to water and oxygen thus causing the community metabolism to switch more often from anoxic to oxic functionalities. Sulfur is essential for building amino acids (e.g., methionine and cysteine), cellular components (sulfolipids) and metabolites (DMSP) (Wasmund et al., 2017). It has a paramount role as electron acceptor (e.g., sulfate and thiosulphate) and electron donor (sulfite and elemental sulfur) (Jørgensen et al., 2019). We found genes related to sulfate-reduction and sulfur-oxidation processes that are coupled in the surface sediment. Overall, sulfur cycling is interconnected with the other elements, the availability of organic matter and redox species present in the lagoon system. In a high-temperature and low-oxygen lagoon scenario, sulfur-related metabolisms might become dominant in the ecosystem functioning.
Phosphorus-cycle-related genes were also highly abundant in our samples, being especially linked to the involvement of nucleic acids in biofilm structures (Karygianni et al., 2020) and organic P degradation. Furthermore, approximately 4.5% of genes (KEGG annotation) in the metagenomes were addressed to peptide degradation. Proteolysis is in fact the fastest hydrolytic activity within surface brackish and coastal sediments, aiming to provide prokaryotes with both organic C and N via amino acids, as highlighted by direct measurements (e.g., Cibic et al., 2012;Franzo et al., 2019). Together with carbohydrates, proteins and lipids are indeed the major constituents of the labile organic matter in the sediments (Pusceddu et al., 1999) and microbes are adapted to utilize these compounds heterotrophically (Oni et al., 2015), thus supporting the concept of sediments as hot spots for organic matter degradation (Lipka et al., 2018).

Cell-Cell Microscale Interactions
Within this framework of reference, we have mined our data set to test hypotheses on the putative functional microscale behavior of the microbes in the lagoon surface sediments (see Figure 9). Starting from DNA, we have asked questions about whether the community was able to display adaptive strategies such as competence and gene exchange to improve their persistence in the environment. Among the mobile genetic elements, what were the most abundant, and was it connected with the site features?
DNA is an integral part of the sediment, carrying genetic information and providing 3D structure by binding to other molecules (Dell'Anno and Corinaldesi, 2004;Corinaldesi et al., 2008;Torti et al., 2015). It is known that viruses are abundant in the sediment and are major structuring agents, keeping the microbial community under check (Danovaro et al., 2001;Breitbart et al., 2004;Breitbart and Rohwer, 2005). Among the microbial defense systems (Makarova et al., 2013), we detected the CRISPR-Cas (Doudna and Charpentier, 2014) associated genes at all the stations, despite the sediment grain size. This suggests that there is a persistent strategy of Bacteria and Archaea to fight off the viral invasions using this adaptive immunity against foreign genetic elements while creating a viral data bank for further viral encounters. Recent studies showed that in Global Ocean Sampling Expedition (GOS) and Tara Ocean data set, CRISPR-Cas were highly diverse and specific to the viral architecture (Cai et al., 2013;Nasko et al., 2019). In addition to CRISPR-Cas, Competence and Sporulation related functions were found throughout the Lagoon. DNA can be taken up during the competence window (Lennon, 2007;Dubnau and Blokesch, 2019), while a population of microbes is undergoing major environmental stress such as oxygen depletion/enrichment or nutrient scarcity. Among that population, a few cells can promote the lysis of the other individuals exploiting the released DNA for ensuring success in the future. Temporally subsequent to competence is the sporulation, an important unique microbial mechanism in marine sediment for surviving burial (Wörmer et al., 2019;Morono et al., 2020).
An ability that microbes use to move information relies on the highly diverse architecture of the mobile genetic elements, the Mobilome (Gillings, 2013;Carr et al., 2021). Prophages and Transposons were lower in Chioggia, possibly due to the site features such as higher water flow and lower TOC that might select against these mobilome strategies. Overall, mapping the sediment mobilome could give insights on rapid adaptation processes to the presence of pollutants in marine habitat due to anthropogenic activities.
Within the Regulation and cell signaling and Nanomachines categories, we have mined and detected genes related to physiology, adaptations and interactions of the microbial communities and consortia in the sediment. Biofilm in the marine environment is a pervasive structuring feature (de Carvalho, 2018). Biofilm formation is essential in keeping the microenvironment stable for the establishing of the microbial communities thus providing protection from grazers and phages, toxic substances (e.g., antibiotics), and accumulating nutrients (Flemming et al., 2016). Toxin-antitoxin (TA) systems regulate dormancy and persistence in microbes (Spoering and Lewis, 2001;Keren et al., 2004;Shah et al., 2006). TA genes are abundant on plasmids, phages, and chromosomes. They are formed by a toxin, that arrests growth by interfering with a vital cellular process, and a cognate antitoxin, that neutralizes the toxin activity to resume normal growth (Page and Peti, 2016). TA plays a major role in supporting microbial resistance to antibiotics by developing a subpopulation of persistent cells. It is interesting to note that in Chioggia, biofilm and TA were lower than the other sites, suggesting that these strategies were not favored in a sandy sediment habitat and or were not so successful as in the other sites. Despite it is known that in sandy sediments, microbes colonize grains via biofilms (Probandt et al., 2017), the biofilm might be more widespread in muddy sediment to create an organic bridge among clay and silt particles.
In terms of competition for space, nutrient, energy sources and ability to resist osmotic stress, phenazine (Price-Whelan et al., 2007;Schiessl et al., 2019), siderophores (Sandy and Butler, 2009;D'Onofrio et al., 2010), two-component systems (Held et al., 2019), osmotic stress genes, and nanomachines for cell-cell warfare were identified in the data set. Microbes constantly under siege could express phenazine, a redox-active pigment produced by Pseudomonas aeruginosa that affect gene expression, metabolic flux, and redox balancing to promote their success. Pseudomonas despite being a multi-drug resistant, biofilm former, quorum sensing human pathogenic microbe is highly versatile and has been found in the marine environment (Kimata et al., 2004). Siderophores are structures produced by microbes that chelate ions, such as Fe to facilitate the acquisition of the insoluble species thus harshening the competition for this element. We identified, within the two-component signal transduction systems (Yamamoto et al., 2005) single transduction for anaerobic respiration, nitrate and nitrite utilization, nitrogen assimilation metals efflux system, fumarate respiration, outer membrane proteins, phosphate assimilation, magnesium transport, and osmotic and envelop stress response, chemotaxis modulation, and flagellar motor switch. This great diversity found in the sediment microbial communities highlighted the adaptive strategies that microbes need to exploit while sensing the external environment at the microscale in order to activate genes that promote their success in the microniche. In the water column, a high abundance of these regulative network systems was found to be linked to copiotrophy and diazotrophy (Held et al., 2019).
The last microscale features we have identified in our dataset were the Type IV (Craig et al., 2019) and Type VI (Cianfanelli et al., 2016) secretion nanomachines, that are involved in cell-cell interactions by allowing nucleic information being shuttle from one individual to the other and attaching with effectors the other organisms. Type IV was found in microbes from acetate-stimulated aquifer sediments (Kantor et al., 2013), thus suggesting conjugation being potentially important to exchange genes and associated metabolic potentials. In the ocean, Type VI secretion systems were found to be globally distributed and associated with particle-attached microbes (Kempnich and Sison-Mangus, 2020).
Finally, motility and chemotaxis patterns, highly abundant in Chioggia, suggested a specific strategy in sandy sediments, that are more permeable to water in comparison to clayey and silty one and where a more diverse microenvironment can develop within the grains.
In sum, microscale-related genes highlight a very elaborated fine structure of interactions and behavior in the surface sediment microbial communities that then results in the large-scale biogeochemical processes on the lagoon system.

Anthropogenic-Related Pressures Influence Benthic Microbial Communities
DNA metabarcoding and metagenomic analyses allowed the exploration of microbial features related to different anthropogenic pressures, of particular importance in the monitoring of fecal and sewage bacteria, antibiotic and metal resistant, hydrocarbon degraders, and xenobiotic degraders microbes.
The use of DNA-based approach is increasingly applied to microbial profiling in wastewater plants, coastal areas, and in tracking and identifying microbial signatures in the environment; such studies contributed to the identification of bacterial taxa associated with sewage-and fecal-contamination (Newton et al., 2013;Tan et al., 2015). Marine sediments are a reservoir of fecal bacteria (Luna et al., 2010) and have the potential to favor the persistence of fecal microbes and to contaminate the overlying water through resuspension, which likely poses important public health and environmental threats (Luna et al., 2016). In our study, microbial signatures showed the highest proportion of sewage-associated taxa in Marghera and the lowest proportion in Palude della Rosa, considered less impacted by human activities. Sacca Sessola and Chioggia had a higher percentage of fecal-associated taxa, being mostly affected by urban contamination (Picone et al., 2016). Overall, we detected several indicators of fecal and sewage contamination (see Figure 6), and their relative abundance and taxonomic assignment are in line with a previous study in central part of the Venice Lagoon (Luna et al., 2016), confirming a diffuse contamination and accumulation of these bacteria in the Lagoon sediment.
Antibiotic resistance spread is one of the main threats to global human, fish, shellfish, and in overall ecosystem health (Yang et al., 2012;Ventola, 2015;Mackenzie and Jeggo, 2019;Helsens et al., 2020). For this reason, the detection, tracking, and assessment of the increasing antibiotic resistant microbes and genes, is necessary in terms of risk management (Heß et al., 2018). The marine environment is considered a global reservoir of antibiotic resistance genes, but we still have little understanding of their presence and diversity in both pelagic and benthic environments (Chen et al., 2013;Hatosy and Martiny, 2015). Several antibiotic resistance genes were detected in our dataset (see Figure 8), with a significantly higher presence in the two more impacted sites of Marghera and Tresse. The Multi-drugs and Beta-lactamase resistance genes were found as the most abundant at all sites, as previously detected in seawater (Hatosy and Martiny, 2015) and sediment (Huang et al., 2019).
Surface sediment metal contamination is another consequence of anthropogenic pressure in transitional environments, representing a potential threat to both public health and marine ecosystems (Han et al., 2011). Metals not only accumulate in the sediment but can be subjected to complex cycling that can make toxic compounds available to benthic organisms; sediment mobilization and resuspension under particular pH and redox conditions may release them in the water column becoming a secondary source of pollution (Turner and Millward, 2002;Eggleton and Thomas, 2004;Zonta et al., 2018). Furthermore, microbes that are antibiotic resistant can carry information for metal resistance and vice versa, thus exacerbating the problem of resistance spreading (Chen et al., 2019). Several toxic compound resistance genes were detected in our dataset, with a significantly higher presence in the two more impacted sites of Marghera and Tresse. These sites presented the highest number of genes related to Cobalt-zinc-cadmium resistance, and to Mercury resistance, due to the release of such metals from industrial plants. In particular, the mercury contamination mainly originates from the chloralkali discharge in Porto Marghera, which occurred until the 1980s (Bloom et al., 2004). Mercury is one of the most severe pollutants, and its flux to Venice Lagoon water is mostly due to the resuspension of contaminated sediment (Bloom et al., 2004), underlying the importance of its regulated management. Chioggia and Palude della Rosa presented the highest number of genes coding for Arsenic resistance. Even if the presence of this contaminant in the lagoon is mainly due to the geochemical characteristics of the sediments, human activities can contribute to its accumulation .
Chioggia presented the highest number of genes related to hydrocarbon degradation. Interestingly, in this site, both surface and sub-surface sediment were recently found to harbor the highest concentrations of PAH and PBC in the whole Venice Lagoon , due to the important urban waste discharge, road/boat traffic, and being located close to the harbor of the second Italian largest fishing fleet.
It was worth noting that among the identified hydrocarbon degraders some taxa were also able to utilize lignin. Lignin is a very complex and hard-to-degrade polymer (Hedges and Mann, 1979) and Bacteria, Archaea and Eukaryotes have evolved adaptive strategies to gain energy from it (Cragg et al., 2015). In the Venice Lagoon, wood has been historically and extensively used for maritime needs, like mooring piles and navigation channel delimitation (named "briccole, " approximately 7,000, formed by two or three poles) (Ghirardini et al., 2010). Briccole are subjected to physical and biotic degradation, the latter mainly due to marine borers (Ghirardini et al., 2010). Since our data highlighted that the sediment represents a reservoir of potential lignin-degrader bacteria, we hypothesize a role of these taxa (alone or in consortia with other microorganisms) in the turnover/recycling of human-introduced wood structures in the Lagoon.

CONCLUSION AND FUTURE PERSPECTIVES
Our analyses showed a clear spatial distribution of the highly biodiverse surface sediment microbial communities in the different sub-basins of the Venice Lagoon, mainly due to low abundant taxa. Salinity, grain size, and TOC concentration were significant drivers in shaping community structure and diversity. Further analysis will be necessary to depict the environmental drivers of microbial dynamics at the small scale.
Metagenomic approach allowed for the first time to assess the potential functions and metabolism of these communities at different levels. At both the ecosystem-and the microscale, subtle but significant differences were detected in the diverse sites, including the main biogeochemical cycles. The focus on anthropogenic-related functions highlighted that chronically polluted areas are hotspots of resistance genes, with implications on sediment management in the Lagoon. Furthermore, our data indicate that an enhanced understanding of microbial microscale features is important to better describe how interactions play out in structuring the community and thus dictating the ecosystem functioning.
In the near future, the metagenome-assembled genomes (MAGs) analysis will give more detail about the link between taxonomy and functions of interest, especially important for the still uncharacterized taxa. This study represents the initial, essential step toward the application of other "omics" approaches in the Venice Lagoon and other transitional sites. Among them, RNA-based approaches (e.g., metatranscriptomics) will benefit from these data to identify the functions and metabolisms that are actually active in the sediment microbial communities. This will allow inferring the microbial response to environmental changes including the increase of temperature, acidification, and extreme flood events.

AUTHOR CONTRIBUTIONS
FM conceptualized the study. PDN acquired the funding. FM and MC performed the management of the project activities.
EB performed the laboratory and bioinformatic analyses, wrote the manuscript with the help, and inputs of all co-authors. All authors contributed to the article and approved the submitted version.

FUNDING
Scientific activity performed with the contribution of the Provveditorato for the Public Works of Veneto, Trentino Alto Adige, and Friuli Venezia Giulia, provided through the concessionary of State Consorzio Venezia Nuova and coordinated by CORILA within Venezia2021 framework. This work benefited from access to SBR-Station Biologique de Roscoff, an EMBRC-FR and EMBRC-ERIC. Financial support was provided by EMBRC-ERIC and the European Marine Research Network (EuroMarine) through their joint 2020 call targeting early career researchers (EB). This work benefited from CINECA-ISCRA C resources (HP10CS2I1S; PI: EB).