Layout 1 INTRODUCTION The search of factors shaping the protistan plankton community has been extended and enhanced in the last years thanks to the combination of molecular data with tools from macroecology. Studies on habitat heterogeneity in different marine ecosystems (Grossmann et al., 2016; Massana et al., 2015; de Vargas et al., 2015; Wu et al., 2020) have shown that plankton microeukaryotes present different community compositions depending both on en- vironmental factors, affecting the habitat features, and on geographical distance on the small scale (Horner-Devine et al., 2004; Martiny et al., 2006; Wu et al., 2020). Habitat heterogeneity is particularly marked in tran- sitional water ecosystems (TWE; European Union, 2000; Elliot and McLusky, 2002), such as coastal lagoons, due to their geomorphology, the close benthic–pelagic coupling, the freshwater inputs, and the connection with the adjacent marine ecosystems through tides and currents, (Basset et al., 2006; McLusky and Elliot, 2007; Vadrucci et al., 2007). These TWE features have a profound effect on plankton composition and distribution at both spatial and temporal scales (Cloern and Jassby, 2010). The Lagoon of Venice (LoV) is the largest Italian la- goon. This TWE belongs to the Long Term Ecosystem Re- search (LTER) national (LTER-Italy), European (LTER- Europe) and global (ILTER) networks, which are essential components of the worldwide efforts to improve our knowledge of the structure and functions of ecosystems and of their long-term response to environmental, societal, and economic drivers (Mirtl et al., 2018). They are distributed networks of research sites for multiple purposes in the fields of ecosystem, critical zone, and socio-ecological research. LTER-Italy consists of 79 research sites, which include ter- restrial, freshwater, transitional and marine ecosystems (Pugnetti et al., 2013; Bergami et al., 2018). For what concerns planktonic protists of the LoV, past LTER activities addressed almost exclusively phyto- plankton communities, for which a consistent knowledge has been gathered through morphology-based studies (Bernardi Aubry et al., 2013, 2017). The LoV community includes taxa that are representative of the complex and open system of the lagoon itself, from the coastal sea to the rivers and the benthic-pelagic system: overall, they show a rather even distribution across distinct lagoon sites, while habitat heterogeneity could barely be detected (Bernardi Aubry et al., 2013). High Throughput Sequencing (HTS) metabarcoding has so far been applied to the LoV bacterial communities (Quero et al., 2017) and, more recently, also to the whole protistan community, in a comparison between ARTICLE Spatial diversity of planktonic protists in the Lagoon of Venice (LTER-Italy) based on 18S rDNA Simona Armeli Minicante,1 Roberta Piredda,2 Stefania Finotto,1 Fabrizio Bernardi Aubry,1* Francesco Acri,1 Alessandra Pugnetti,1 Adriana Zingone2 1Institute of Marine Sciences, National Research Council (ISMAR-CNR), Venice; 2Department of Integrative Marine Ecology, Stazione Zoologica Anton Dohrn, Naples, Italy ABSTRACT Transitional waters are subject to a high degree of variability in space and time. In this study, protist plankton communities of the Lagoon of Venice were compared among four sites characterised by different environmental conditions with a metabarcoding approach. High throughput sequencing (HTS) of the V4-18S rDNA fragment in 32 samples collected on four dates, from April 2016 to February 2017, produced 1,137,113 reads, which were grouped into 4058 OTUs at 97% similarity. Bacillariophyta and Ciliophora were the most abundant groups in the entire dataset in terms of read number (27.6% and 16.6%, respectively), followed by Dinophyta (10.9%), Cryptophyceae (9.7%), and Syndiniales (6.1%). The contribution of protist groups markedly varied across the seasons, but spatial differences were also recorded in the lagoon. In April, a higher contribution of Bacillariophyta characterized St1 and 5 (68.0% and 61.1%), whereas Sts2 and 3 showed a higher percentage of Ciliophora (18.6 and 23.4%, respectively) and dinoflagellates (10.3 and 7.7%). In July, diatom blooms occurred at Sts1, 2 and 3, with some differences in the dominant species. At St2 Dinophyta reached the highest contribution of the whole sampling period in the area (30.6%), while St5 was quite distinct, with a low contribution of diatoms and a dominance of Ciliophora (34.0%) and Trebouxiophyceae (36.4%). All the stations in No- vember were characterized by relatively high abundance of Ciliophora (21.4-51.9%). In February, diatom contribution was relevant only at St5 (29.3%), Teleaulax acuta peaked at St3 (ca. 36%), Syndiniales at St2 (38.8%) and Dictyochophyceae at St1 (24.2%). The α-diversity indexes (observed OTUs, Shannon and Pielou evenness) showed a high variability over space and time. Diversity and community composition were rather similar between the intermediate and deeper Sts2 and 3 on all sampling dates whereas they at time differed between the landward and shallow Sts1 and 5. While the most marked differences occurred over the temporal scale, the depth of the station and the relatedness with the external marine coastal environment appear to play a major role in the spatial distribution of protist communities within the lagoon. Overall, the quite stable spatial differences in spite of the significant seasonal changes reflected the hydro-geological heterogeneity of the sampling stations, indicating a major influence of the land- sea gradient in the lagoon. No n- co mm er cia l u se on ly S. Armeli Minicante et al.36 the lagoon and the adjacent marine waters of the Gulf of Venice (Armeli Minicante et al., 2019). In spite of potential pitfalls of the metabarcoding approach, which can poorly resolve the diversity of some taxa (Piredda et al. 2018) or miss them altogether (Massana et al., 2015), the latter study based on the V4-18S rRNA metabarcodes largely increased the diversity knowledge not only for protists that have tra- ditionally been neglected (i.e., heterotrophs, parasites, pi- coeukaryotes and other featureless groups), but also for the main phytoplankton taxa studied in the long term with mor- phology-based approaches (i.e., diatoms and dinoflagel- lates). In addition, HTS results highlighted profound differences in the structure of the protistan communities be- tween the lagoon and the external coastal waters of the Gulf of Venice: heterogeneity appears strong enough to allow for ecological segregation in the two environments, despite no clear barrier to dispersal processes among local protist communities. In this study, we explore the effects of habitat het- erogeneity within the LoV by comparing protistan com- munities among four individual LoV sites that were lumped together in the previous HTS study. The main goal of this work is to test whether the overall picture of rather spatially homogeneous LoV communities emerging from morphologically-based studies on phytoplankton is con- firmed by the higher resolution and better coverage al- lowed by the HTS-metabarcoding approach. We aim at assessing if environmental selection plays a more impor- tant role than dispersal in shaping the microeukaryotic communities within the LoV and identifying the most rel- evant local conditions that could affect them. METHODS Study area The LoV (Northern Adriatic Sea, Mediterranean Sea) is a microtidal, polyhaline lagoon classified as a tran- sitional water body (European Union, 2000) and is the largest wetland (550 km2) in the Mediterranean Sea (Pog- gioli, 2008). This delicate and fragile ecosystem is sur- rounded by densely inhabited and industrial areas and is affected by high numbers of tourists, as well as intensive fisheries and aquaculture activities. It has an average depth of 1 m and is morphologically characterized by the presence of large shallow areas and a network of deeper (5–10 m) channels. Three inlets interrupt the sandbars separating the lagoon from the sea, allowing water turnover during tidal cycles, whose amplitude is 100 cm with maxima of 150 cm. Water residence time, resulting from the interactions of tide, wind and topography, ranges from a few days close to the inlets to one month in land- ward areas (Umgiesser et al., 2014; Ghezzo et al., 2015). Twelve main tributaries discharge an annual average of about 35 m3 s−1 of freshwater into the lagoon, with sea- sonal peaks in spring and autumn. In recent years aquatic angiosperms have recolonized the lagoon, especially in the central and northern basins (Sfriso and Buosi, 2018); however, in large areas the main primary producer is phy- toplankton (Acri et al., 2004, Sfriso et al., 2005; Bernardi Aubry et al., 2013). In this study four sampling campaigns were con- ducted in April, July, November 2016 and February 2017, respectively, at four stations (St1, St2, St3 and St5, Fig.1) that are regularly sampled for phytoplankton and abiotic factors within LTER activities (Bernardi Aubry et al., 2013). The stations are all influenced, at different extents, by marine and freshwater inputs: St1 (S. Giuliano, 2 m depth ) is in an area of intense maritime traffic (Bianchi et al. 1996, Socal et al. 1999) which collects urban waste from the town of Mestre; St2 (Marghera, 9.5 m) is af- fected by industrial pollution (Perin 1975, Guerzoni et al. 2007); St3 (Fusina, 3.5 m) is influenced by heat emissions from the Porto Marghera power station (Alberighi et al., 1992, Socal et al., 1999); St5 (Palude della Rosa, depth 2.7 m) is a typical inland marshy lagoon area (Bianchi et al., 1999). Differences among the stations mainly depend on the varying degrees of impact of the interactions among the adjacent marine waters, riverine inputs and benthic-pelagic coupling (Bernardi Aubry et al., 2013). Abiotic parameters For each station, temperature and salinity were measured with a CTD SBE 911 probe and water samples were collected at the near-surface layer using a Niskin bottle. For chlorophyll a (Chl a), 500 mL of water samples were immediately filtered through Whatman GF/F fibre- glass filters (nominal porosity = 0.7 µm), which were stored frozen and subsequently analyzed according to Holm-Hansen et al. (1965). Dissolved Inorganic Nitrogen (DIN), as sum of ammonium (N-NH4), nitrites (N-NO2) and nitrates (N-NO3), was analyzed along with orthophos- phates (P-PO4) and orthosilicates (Si-SiO4) with a Systea EasyChem Plus according to Grasshoff et al. (1983). Filtration, DNA extraction and sequencing At each station, 3 L of seawater were prefiltered on a 200 μm mesh-size net and then filtered onto cellulose ester 1.2 μm pore size filters (47 mm Ø, Whatman) using a peristaltic pump. For each station a duplicate was per- formed, obtaining 32 filters that were stored at -80°C until molecular analysis. Total DNA from each filter was ex- tracted using the DNeasy 96 Plant Kit (QIAGEN) and concentrations were determined with the Qubit dsDNA HS kit (Thermofisher). DNA samples were stored at −80°C until PCR. The hypervariable V4 region of eukary- No n- co mm er cia l u se on ly Planktonic protist diversity in the Lagoon of Venice 37 ote SSU rDNA gene was amplified and sequenced (2 X 250 bp sequencing) on an Illumina MiSeq platform as de- scribed in Piredda et al. (2017). Sequence analyses Paired-end reads were processed using Mothur v.1.33.0 (Schloss et al., 2009). Contigs between read pairs were assembled and differences in base calls in the overlapping region were solved using ΔQ parameter as described in Kozich et al. (2013). Primer sequences were removed (pdiffs=3) and no ambiguous bases were al- lowed; the maximum homopolymer size was 8 bp. The remaining sequences were de-replicated and screened for chimeras using UCHIME in de novo mode (Edgar et al., 2011). Sequences were clustered into operational taxonomic units (OTUs) at 97% of similarity using vsearch (Rognes et al., 2016) clustering (method=dgc) through mothur. OTUs containing only one read (single- ton) were removed from downstream analyses. Taxo- nomic assignment was performed on a single representative sequence from each OTU (the most abun- dant) using BLASTN (Altschul et al., 1990) against the PR2 database (v.4.10.0; https://doi.org/10.6084/m9. figshare.5913181; Guillou et al., 2013), discarding the assignations with similarity ≤90% and query coverage ≤70% of the sequence length. Statistical analyses The vegan R package (R Core Team, 2014; Oksa- nen et al., 2016) was used for multivariate analyses. The dataset was normalized with a random subsam- pling to the second-lowest number of sequences (n = 74,396) with the ‘rrarefy’ function. Observed OTUs and α-diversity estimators (Shan- non diversity and Pielou evenness) were calculated and Non-metric multidimensional scaling (NMDS) was per- formed using the ‘metaMDS’ function (vegan) based on a Bray-Curtis dissimilarity matrix. PerMANOVA (adonis function in vegan) was used to detect significant associa- tion between environmental parameters and protist com- munity structure. The selected variables were used to perform Canonical Correspondence Analysis (CCA). RESULTS Abiotic parameters and chlorophyll a Abiotic parameters and Chl a showed a wide vari- ability over time and among stations (Tab. 1). Salinity val- ues (21-32.3) most of the times were lowest at Sts1 and 5 (>21 and 23, respectively), while they were above 30 all year round at the other two stations. Temperature ranged between 9.4°C (St2, February) and 31.5°C (St3, July), with Fig. 1. Study area with the four sampling stations. No n- co mm er cia l u se on ly S. Armeli Minicante et al.38 marked differences (6-10°C) among stations on the same sampling date, and St3 always showing the highest values. DIN and silicate concentrations peaked at all the stations, and particularly at Sts1 and 5, in November (25.4-88.1 µM) and February (9.1-61.7 µM), the maxima occurred mostly at St1 in every date. Phosphates were always highest at St1 (0.6-3 µM), while at the other stations they fluctuated be- tween 0,1 and 1.1 µM without a clear seasonal or spatial pattern. Chl a concentrations were higher at Sts1 and 5 (2.7 and 2.2 µgL–1, respectively) in April, then they peaked in July at Sts1, 2 and 5 (12.6, 11.8 and 9.6 µgL–1, respectively) and were relatively low at all stations in November (0.2- 0.7 µgl–1) and February (0.1-0.7 µgL–1). The highest Chl a concentrations were recorded at St1 on all dates but Febru- ary, when the maximum was at St5. Protist diversity and seasonality The total cleaned dataset contained 3,588,416 reads from 32 samples. Duplicate samples were generally cou- pled in hierarchical clustering based on Bray-Curtis dis- tance (data not shown) and were hence pooled into 16 samples. After the normalization, the dataset consisted of 1,137,113 V4-18S reads (Supplementary Material 1), clustering into 4,058 OTUs mostly belonging to 34 high- level taxonomic groups (Fig. 2). Bacillariophyta and Ciliophora were the most abundant groups in terms of read number (27.6% and 16.6%, respectively), followed by Dinophyta (10.9%), Cryptophyceae (9.7%), Syndiniales (6.1%). Chryso- phyceae, Dictyochophyceae, Mamiellophyceae, Tre- bouxiophyceae, MAST (nanoheterotrophic marine stramenopiles) and Chlorophyceae ranged between 2 and 5%. The remaining 23 taxonomic groups were included in “other Eukaryotes” (8.9%). Over the four sampling dates, a change in the pro- portion of the main groups was evident at all stations, which at times also differed among them (Fig. 3). Common to Sts1, 2 and 3 was the high relative abundance of Bacillar- iophyta both in April (max 67.9%, at St1) and July (max 87%, at St3) and their very low abundance in November (1.1%, at St2) and February (1.9%, at St2). At St5, Bacil- Tab. 1. Environmental variables at the four stations of the Lagoon of Venice investigated on four sampling dates. Month Station Temperature Salinity DIN Si-SiO4 P-PO4 Chl-a (°C) (µM) (µM) (µM) (µgL–1) April ST1 14.1 20.9 30.3 36.7 0.5 2.7 ST2 17.3 30.1 17.9 13.1 0.5 0.8 ST3 24.3 29.5 23.1 18.3 0.3 0.1 ST5 14.1 30.3 13.9 10.8 0.1 2.2 July ST1 25.3 29.5 15.2 41.6 2.3 12.5 ST2 28.8 31.3 6.5 13.6 0.2 11.7 ST3 31.5 30.0 17.5 16.6 1.2 1.7 ST5 25.5 27.1 12.3 43.8 0.6 9.5 November ST1 14.5 28.3 83.9 88.1 2.5 0.7 ST2 14.7 31.0 36.9 28.4 1.0 0.3 ST3 22.1 31.6 36.0 25.3 0.8 0.2 ST5 14.0 25.9 63.8 37.2 0.4 0.5 February ST1 10.4 21.0 75.8 61.6 2.9 0.5 ST2 9.4 31.3 53.8 27.1 0.7 0.1 ST3 15.6 32.9 27.4 12.7 0.4 0.5 ST5 11.8 23.3 24.4 9.1 0.2 0.7 Fig. 2. Taxonomic diversity of the normalized dataset as percent- age of reads belonging to the main 10 high-level taxonomic groups and to the broad group “Other Eukaryotes” which in- cludes: Apusozoa, Radiozoa, Amoebozoa, Stramenopiles, Api- complexa, Perkinsea, Palagophyceae, Bolidophyceae, Telonemia, Pirsonia, Bicoecea, MOCH, Prymnesiophyceae, Centroheliozoa, Katablepharidaceae, Choanoflagellatea, Labyrinthulea, Prasino- phyceae, Oomycota, Chlorodendrophyceae, Chlorophyceae, Cer- cozoa, Rhodophyta, Picozoa. No n- co mm er cia l u se on ly Planktonic protist diversity in the Lagoon of Venice 39 lariophyta contribution was highest in April (61.1%) too, but low in July (10.9%) and relatively high again in Febru- ary (29.3%). Other seasonal patterns shared among stations concerned the higher contribution of Trebouxiophyceae in July (up to 23.6% and 5 36.4% at Sts1 and 5, respectively), of Ciliophora at all stations in November (21.3-51.9%), and the low contribution of Cryptophyceae and Mamiel- lophyceae in July. Besides, at all stations Dictyochophyceae were more abundant at almost all stations in February (<24.1% at St1), while the contribution of Chrysophyceae was significant in November, especially at St5 (12.9%). Dinophyta attained more than the 20% of the total abun- dance only at St2 in July (30.6%) and November (23.7%). Within each sampling date, some dissimilarities among stations were also evident. In April, the main dif- ference was the highest contribution of Bacillariophyta at Sts1 and 5, which also showed the highest Chl a concen- trations. Yet, diatom species composition differed between these two stations: Thalassiosira concaviuscula domi- nated both the group and the whole protist assemblage at St1 (56.6% of total reads), while it covered less than 1% at St5, where the dominant species (50%) was Cyclotella sp. At the other two stations, diatoms were less abundant and mainly represented by Thalassiosira spp., Chaeto- ceros tenuissimus, Navicula and Skeletonema. Sts2 and 3 differed from Sts1 and 5 also for the highest abundance of Ciliophora (18.6% and 23.4%), mainly represented by Strombidida and by the tintinnid Tintinnopsis, and of di- noflagellates (<10%), with the autotrophic species Gymnodinium dorsalisulcum and Heterocapsa pygmaea. Mamiellophyceae, mainly Ostreococcus mediterraneus, were most abundant at St3 (10.9%), Cryptophyceae at St2 (15.0%) with Teleaulax acuta and Rhodomonas sp. In July a diatom bloom occurred at Sts1, 2 and 3, with the highest relative percentage contribution (87.5 %) at St. 3 and the lowest (50.7%) at St1. Chaetoceros tenuis- simus dominated at Sts1 and 3, the cymatosiracean Min- utocellus polymorphus at St2. There, Dinophyta, with G. dorsalisulcum and H. pygmaea (92.0%), attained the highest contribution (30.6%) of the whole sampling pe- riod. St5 was quite distinct from the others due to a much lower percentage of diatoms (10.9%) and the dominance of Ciliophora (34%, mainly Parastrombidinopsis sp.) and Trebouxiophyceae (36.4%), the latter mainly represented by the green coccoid Picochlorum which was also abun- dant at St1 (23.6%). In November all the stations were characterized by the high abundance and diversity of Ciliophora (21.3- 51.9%), among which the most abundant were an un- known Strobilidiidae, Rimostrombidium and Strombidium. In addition to Ciliophora, Sts2 and 3 showed a quite similar composition in terms of both groups and species, with a high proportion of Dinophyta (23.7% at St2, 15.1% at St3), mainly represented by H. pygmaea, the parasitic group Syndiniales (13.3 at St2, 16.1% at St3), heterotrophs, such as MAST (6.7% at St2, 4.4% at St3) and other eukaryotes (such as Picozoa and other non-photosynthetic stramenopiles, 8.9% and 9.3% Fig. 3. Composition of the protist community at the four stations at each sampling date, based on the normalized dataset. No n- co mm er cia l u se on ly S. Armeli Minicante et al.40 at Sts 2 and 3, respectively ). Cryptophyceae (with a high proportion of Teleaulax acuta) peaked at St5 (20.5%) to- gether with the Dictyochophyceae Pedinellales (18.7%). In February St5 differed from all the others for the dominance of diatoms (29.3%) mainly represented by ben- thic taxa, either marine (Undatella and Achnantes) or fresh- water (Cymbella), a higher percentage of ciliates (15.0%) and the unicellular Rhodophyta (Rhodella sp.) as the most abundant taxon (10.9%). The dominance of Teleaulax acuta was remarkable at St3 (36%), while Syndiniales were abun- dant at St2 (38.8 %) and Dictyochophyceae at St1 (24.2%, mainly Apedinella radians). Dinophyceae at St2 and St3 were mainly represented by H. pygmaea. Trends of alpha diversity (Fig. 4) showed a great variability at both time and space scales. The effect of the diatom blooms was more evident in July, resulting in the lowest diversity values particularly at Sts 2 and 3. In Feb- ruary, observed OTU number was the lowest at St5, which however showed the highest evenness value. Overall, di- versity indexes were often similar within the two couples of stations, St1-St5 and St2-St3. The NMDS ordination highlighted a relatively reg- ular pattern in relation to the spatial and temporal scales (Fig. 5). NMDS2 values below or above zero split Sts2 and 3 (bottom of the plot) from Sts1 and 5 (top), whereas NMDS1 separated February and November communities (to the left) from April and July ones (to the right), the lat- ter being the most peculiar ones at all four sites. Sts2 and 3 were strictly coupled in all sampling events but July, whereas Sts1 and 5 were slightly more distant but still coupled in all sampling but February. PerMANOVA selected three environmental vari- ables significantly correlated with the community com- position: DIN, salinity, and temperature. In the CCA (Fig. 6) performed using these variables, the first two canonical axes only explained 23% of the total variance. Yet, the analysis confirmed the clear separation of the July com- munity, which showed a high correlation with tempera- ture. February and November communities were associated with high DIN values. Sts2 and 3 were coupled in all seasons, although at times more distant than in the NMDS, whereas Sts1 or 5 were often closer to one of the other two stations than between them. DISCUSSION High spatial and temporal variability is the hall- mark for coastal environments and especially in TWEs, which are sensitive to an array of external drivers and pressures. In particular, in the LoV the low average depth causes an accelerated exchange and biogeochemical cy- cling between the sediments and the water column (Soli- doro et al., 2010), which interacts with the variable Fig. 4. Trend of α-diversity indexes: observed OTUs, Shannon diversity and Pielou evenness. No n- co mm er cia l u se on ly Planktonic protist diversity in the Lagoon of Venice 41 nutrient load from land, tide cycle and meteo-climatic-dri- ven hydrodynamism producing a complex space-time dy- namic. To date the extent to which this complexity is reflected in the microbial community of the lagoon has only partially been understood also because of the limited resolution and partial coverage of the microbial diversity permitted by classical, microscopy-based approaches. In this study, the V4 -18S rRNA metabarcode analysis al- lowed a deeper and more complete assessment of mi- croeukariote communities of the area. Although based on only four time points, the con- ditions recorded in this study can be considered quite in- dicative of the four seasons as they emerge from the previous knowledge of the area. The low plankton abun- dance during the autumn-winter time, the vegetative re- covery in winter-spring and the summer peak match the “typical” seasonal trend observed for phytoplankton com- munities during 10 years of observations in the LoV (Bernardi Aubry et al., 2013). The dominance of Bacillar- iophyta in summer (and, of lesser importance, in spring) and their low contribution in autumn-winter, when our sampling may have missed the typical Skeletonema peak, are also confirmed by light microscopy analysis (data not shown). While the dominance of tiny Chaetoceros in spring–summer period was already known in the area (Bernardi Aubry et al., 2013), the importance of the cy- matosyracean Minutocellus polymorphus at St3 was not recognized before. The highest contribution of Dinophyta in summer was already reported in literature but with lower abundance (average 2%; Bernardi Aubry et al., 2013). However, dinoflagellate abundance can be overes- timated by HTS studies because of the large amount of multicopy genes (Piredda et al., 2017; Gong et al., 2019). Several species commonly found as abundant or featuring different seasons in this study were also reported as char- acterizing in others lagoons or estuarine environments, Fig. 6. Canonical Correspondence Analysis (CCA) based on the normalized dataset, with the three environmental variables selected through PerMANOVA. Fig. 5. Non-metric MDS ordination of protist community com- position at the four stations, based on the Bray-Curtis dissimi- larity matrix calculated on the OTU table. Stress value = 0.099. No n- co mm er cia l u se on ly S. Armeli Minicante et al.42 such as Teleaulax acuta, Cyclotella sp. and Minutocellus polymorphus (Sarno et al., 1993; Bérard-Therriault et al., 1999; Melo et al., 2010; Guiry and Guiry, 2019). From the diversity standpoint, spatial and temporal variations of het- erotrophic, mixotrophic, parasitic such groups as Cilio- phora, MAST and Syndiniales, as well as of autotrophic and heterotrophic picoeukaryotes (e.g., Trebouxiophyceae and Picozoa) emerging in this study represent an absolute novelty and hint at quite intricate and dynamic trophic in- teractions within the microbial food web of the lagoon. For what concerns spatial differences, they appeared quite evident for the whole microeukaryotic communities compared to the overall homogeneity emerging from pre- vious studies that only addressed phytoplankton. In those studies, despite some weak differences among areas, the most abundant species were shared by all stations investi- gated, whereby community composition appeared to be prevalently shaped by the general ecological characteristics of the lagoon ecosystem itself rather than by the specific characteristics of individual areas (Acri et al., 2004; Bernardi Aubry et al., 2004, 2013). The latter picture is based on a much higher number of observations, but the higher resolution allowed by metabarcoding may depict more effectively differences and commonalities among sites. Interestingly, some regularities in the spatial hetero- geneity persisted across the sampling dates, i.e., in the cou- pling between Sts1 and 5 and between Sts2 and 3, where the structure of the protist communities probably reflects different hydro-geological characteristics. Sts2 and 3 are characterized by greater depths (11 m and 4 m, respec- tively) compared to Sts1 and 5 (1.5 m and 2.5 m) which, being the most landward sites, are more affected by resus- pended (e.g., pennate diatoms) and/or freshwaters species from the mainland than by the incoming tide from the in- lets. The rather regular coupling between Sts2 and 3 may be explained by their position in the nearness of the oil tanker channel (Canale dei Petroli) that, connecting with the Malamocco inlet, is more influenced by marine waters from the Gulf of Venice. This area is characterised by higher transit time and lower residence time respect to the other two sampling stations (Guerzoni and Tagliapietra, 2006). Instead, Sts1 and 5 are characterized by a low seabed and by a higher contribution of resuspended species (e.g., diatoms), and can be affected to varying degrees by the in- fluence of the inland anthropogenic inputs (e.g., terrestrial waste in St1 or fresh water in St5), which are also very vari- able over time. This latter interpretation is confirmed by the higher abundance, at Sts2 and 3, of taxa typical of the ex- ternal marine environment such as Cryptophyceae, Dino- phyta, Syndiniales and MAST (Armeli Minicante et al., 2019) with respect to Sts1 and 5. The overall protistan community variance could be related to the average environmental conditions of the sampling stations, in turn depending on the site location in the lagoon, rather than to the instantaneous measured values of the environmental variables during each sam- pling event, which actually explained only a minor part of the variance (23%). In conclusion, this HTS metabarcoding investiga- tion reveals relatively stable spatial differences and sub- stantial seasonal changes in the structure of the whole protistan communities of the LoV, where habitat hetero- geneity and connection coexist. Protistan groups that have distinct functional roles (e.g., Ciliophora, MAST, Syn- diniales and picoeukaryotes), so far never investigated in the LoV, also showed marked differences over the tem- poral and spatial scale, hinting at complex and dynamic relationships within the microbial community which war- rant further investigations. Overall some different hydro-geological character- istics that typify different areas contribute consistently in shaping the spatial microeukaryotic community configura- tion. In particular, the depth of the stations, which modu- lates the benthic-pelagic coupling, and the degree of relatedness with the external marine coastal environment or with the mainland, appear to play a major role in the spa- tial distribution of protist community within the lagoon. ACKNOWLEDGMENTS This work was possible thanks to partial financial support by Italy-LifeWatch-ERIC. RP was supported by the Italian MIUR Flagship project RITMARE and the Eu- ropean Union’s Horizon 2020 Research and Innovation Program EMBRIC (GA 654008). Corresponding author: fabrizio.bernardi@ismar.cnr.it Key words: 18S rDNA; Lagoon of Venice; transitional environ- ments; protists plankton diversity; habitat heterogeneity; LTER- Italy. Received: 21 May 2020. Accepted: 23 June 2020. This work is licensed under a Creative Commons Attribution Non- Commercial 4.0 License (CC BY-NC 4.0). ©Copyright: the Author(s), 2020 Licensee PAGEPress, Italy Advances in Oceanography and Limnology, 2020; 11:8961 DOI: 10.4081/aiol.2020.8961 REFERENCES Acri F, Bernardi Aubry F, Bianchi F, Boldrin A, Comaschi A, Rabitti S, Socal G, 2004. Changes in nutrients and plankton communities in the lagoon of Venice from 1970s until today. J. Mar. Sys. 51: 321-329. No n- co mm er cia l u se on ly Planktonic protist diversity in the Lagoon of Venice 43 Alberighi L, Bianchi F, Cioce F, Socal G, 1992. [Osservazioni du- rante un bloom di Skeletonema costatum in prossimità della centrale termoelettrica ENEL di Fusina Porto-Marghera (Venezia)].[Article in Italian]. Oebalia 17:321-322. Armeli Minicante S, Piredda R, Quero GM, Finotto S, Bernardi Aubry F, Bastianini M, Pugnetti A, Zingone A, 2019. Habitat heterogeneity and connectivity: effects on the planktonic protist community structure at two adjacent coastal sites (the Lagoon and the Gulf of Venice, Northern Adriatic Sea, Italy) revealed by metabarcoding. Front. Microbiol. 10:2736. doi: 10.3389/fmicb.2019.02736. Basset A, Sabetta L, Fonnesu A, Mouillot D, Do Chi T, Viaroli P, Giordani G, Reizopoulou S, Abbiati M, Carrada GC, 2006. Typology in Mediterranean transitional waters: new challenges and perspectives. Aquatic Conserv: Mar. Freshw. Ecosyst. 16:441-455. doi: 10.1002/aqc.767. Bergami C, L’Astorina A, Pugnetti P, 2018. [I cammini della rete LTER-Italia - Il racconto dell’ecologia in cammino].[Book in Italian]. CNR Edizioni, Rome; 180 pp. Bernardi Aubry F, Berton A, Bastianini M, Socal G, Acri F, 2004. Phytoplankton succession in a coastal area of the NW Adriatic over a 10-years sampling period (1990-1999). Con- tinent. Shelf Res. 24:97-115. Bernardi Aubry F, Acri F, Bianchi F, Pugnetti A, 2013. Looking for patterns in the phytoplankton community of the Mediter- ranean microtidal Venice Lagoon: evidence from ten years of observations. Sci. Mar. 7:47-70. doi: 10.3989/scimar. 03638.21. Bernardi Aubry F, Pugnetti A, Roselli L, Stanca E, Acri F, Finotto S, Basset A, 2017. Phytoplankton morphological traits in a nutrient-enriched, turbulent Mediterranean mi- crotidal lagoon. J. Plankton Res. 39:564–576. doi: 10.1093/plankt/fbx008. Cloern JE, Jassby AD, 2010. Patterns and Scales of Phytoplank- ton Variability in Estuarine–Coastal Ecosystems. Estuar. Coast. 33: 230. https://doi.org/10.1007/s12237-009-9195-3. de Vargas C, Audic S, Henry N, Decelle J, Mahé F, Logares R., Lara E, Berney C, Le Bescot N, Probert I, Carmichaell M, Poulain J, Romac S, Colin S, Aury JM, Bittner L, Chaffron S, Dunthorn M, Engelen S, Flegontova O, Guidi L, Horák A, Jaillon O, Lima-Mendez G, Lukeš J, Malviya S, Morard R, Mulot M, Scalco E, Siano R, Vincent F, Zingone A, Dimier C, Picheral M, Searson S, Kandels-Lewis S, Tara Oceans Coordinators, Acinas SG, Bork P, Bowler C, Gorsky G, Grimsley N, Hingamp P, Iudicone D, Not F, Ogata H, Pe- sant S, Raes J, Sieracki ME, Speich S, Stemmann L, Suna- gawa S, Weissenbach J, Wincker P, Karsenti E, 2015. Eukaryotic plankton diversity in the sunlit ocean. Science 348:1261605. doi: 10.1126/science.1261605. Edgar RC, Haas BJ, Clemente JC, Quince C, Knight R. 2011, UCHIME improves sensitivity and speed of chimera detec- tion. Bioinformatics 27:2194-2200. doi: 10.1093/bioinfor- matics/btr381 Elliott M, McLusky DS, 2002. The need for definitions in un- derstanding estuaries. Estuar. Coast. Shelf Sci. 55: 815-827. doi: 10.1006/ecss.2002.1031. European Union, 2000. Common Implementation strategy for the Water Framework Directive (2000/60/EC) - Guidance Document No 5 Transitional and Coastal Waters – Typology, Reference Conditions and Classification Systems. Office for Official Publications of the European Communities; Lux- embourg. Ghezzo M, De Pascalis F, Umgiesser G, Zemlys P, Sigovini M, Pérez-Ruzafa CMA, 2015. Connectivity in three European coastal lagoons. Estuaries Coast. 38:1764-1781. doi: 10.1007/s12237-014-9908-0. Gong W, Marchetti A, 2019. Estimation of 18S gene copy num- ber in marine eukaryotic plankton using a next-generation sequencing approach. Front. Mar. Sci. 6. doi: 10.3389/fmars. 2019.00219 Grasshoff K, Ehrhardt M., Kremling K, 1983. Methods of sea- water analysis. Wiley-VCH Verlag; pp. 419. doi: 10.1002/iroh.19850700232.. Grossmann L, Jensen M, Heider D, Jost S, Glücksman E, Har- tikainen H, Mahamdallie SS, Gardner M, Hoffmann D, Bass D, Boenigk J, 2016. Protistan community analysis: key find- ings of a large-scale molecular sampling. ISME J. 10:2269- 2279. doi: 10.1038/ismej.2016.10. Guerzoni S, Tagliapietra D. 2006. [Atlante della laguna. Venezia tra terra e mare].[Book in Italian]. Marsilio, Venice: 242 pp. Guerzoni S, Rossigni P, Sarretta A, Raccanelli S, Ferrari G, Molinaroli E, 2007. POPs in the Lagoon of Venice: budgets and pathways. Chemosphere 67:1776-1785. Holm-Hansen O, Lorenzen CJ, Holmes RW, Strickland JDH, 1965. Fluorometric determination of chlorophyll. J. Conseil Perm. Int. Explor. Mer. 30:3-15. Horner-Devine M, Lage M, Hughes J, Bohannan BJM, 2004. A taxa–area relationship for bacteria. Nature 432:750–753. doi: 10.1038/nature03073. Kozich JJ, Westcott SL, Baxter NT, Highlander SK, Schloss PD, 2013. Development of a dual-index sequencing strategy and curation pipeline for analyzing amplicon sequence data on the MiSeq Illumina sequencing platform. Appl. Environ. Microbiol. 79:5112-5120. doi: 10.1128/AEM.01043-13 Martiny J, Bohannan B, Brown J, Colwell RK, Fuhrman JA, Green JL, Horner-Devine MC, Kane M, Krumins JA, Kuske CR, Morin PJ, Naeem S, Ovreås L, Reysenbach AL, Smith VH, Staley JT, 2006. Microbial biogeography: putting mi- croorganisms on the map. Nat. Rev. Microbiol. 4:102-112. doi: 10.1038/nrmicro1341. Massana R, Gobet A, Audic S, Bass D, Bittner L, Boutte C, Chambouvet A, Christen R, Claverie JM, Decelle J, Dolan JR, Dunthorn M, Edvardsen B, Forn I, Forster D, Guillou L, Jaillon O, Kooistra WH, Logares R, Mahé F, Not F, Ogata H, Pawlowski J, Pernice MC, Probert I, Romac S, Richards T, Santini S, Shalchian-Tabrizi K, Siano R, Simon N, Stoeck T, Vaulot D, Zingone A, de Vargas C, 2015. Ma- rine protist diversity in European coastal waters and sedi- ments as revealed by high�throughput sequencing. Environ. Microbiol. 17:4035-4049. doi: 10.1111/1462-2920.12955 McLusky DS, Elliott M, 2007. Transitional waters: A new ap- proach, semantics or just muddying the waters? Estuar. Coast. Shelf Sci. 71:359-363. doi: 10.1016/j.ecss.2006.08.025 Mirtl M, Borer E, Djukic I, M.Forsiusc, Haubolda H, Hugod W, Jourdane J, Lindenmayerf D, McDowellg WH, Muraokah H, Orensteini DE, Pauwd JC, Peterseila J, Shibataj H, Wohnera C, Yuk X,Haaseel P, 2018. Genesis, goals and achievements of long�term ecological research at the global scale: a critical review of ILTER and future directions. Sci. Tot. Environ. 626:1439– 62. No n- co mm er cia l u se on ly S. Armeli Minicante et al.44 Oksanen J, Blanchet FG, Kindt R, Legendre P, Minchin PR, O’Hara RB, Simpson GL, Solymos P, Stevens MHH, Szoecs E, Wagner H, 2016. vegan: community ecology package, version 2.4-1. Available from: https://cran. r-project. org/web/packages/vegan/index. html Percopo I, Ruggiero MV, Balzano S, Gourvil P, Lundholm N, Siano R, Tammilehto A, Vaulot D, Sarno D, 2016. Pseudo-nitzschia arctica sp. nov., a new cold-water cryptic Pseudo-nitzschia species within the P. pseudodelicatissima complex. J. Phycol. 52:184-99. doi: 10.1111/jpy.12395 Perin G, 1975. [L’inquinamento chimico della Laguna di Venezia, p. 47-89]. In: Cons. Depur. Acque della Z. I. Porto Marghera (ed.), [Problemi dell’inquinamento lagunare].[in Italian]. Venezia. Piredda R, Tomasino MP, D’Erchia AM, Manzari C, Pesole G, Montresor M, Kooistra WH, Sarno D, Zingone A, 2017. Diversity and temporal patterns of planktonic protist assem- blages at a Mediterranean LTER site. FEMS Microbiol. Ecol. 93. doi: 10.1093/femsec/fiw200 Piredda R, Claverie JM, Decelle J, deVargas C, Dunthorn M, Edvardsen B, Eikrem W, Forster D, Kooistra WHCF, Loga- res R, Massana R, Montresor M, Not F, Ogata H, Pawlowski J, Romac S, Sarno D, Stoeck T, Zingone A, 2018. Diatom diversity through HTS-metabarcoding in coastal European seas. Sci. Rep. 8:18059. doi: 10.1038/s41598-018-36345-9 Power ME, Mills LS, 1995. The keystone cops meet in Hilo. Trends Ecol. Evol. 10:182-184. doi: 10.1016/S0169- 5347(00)89047-3 Pugnetti A, Acri F, Bernardi Aubry F, Camatti E, Cecere E, Facca C, Franzoi P, Keppel E, Lugliè A, Mistri M, Munari C, Padedda BM, Petrocelli A, Pranovi F, Pulina S, Satta CT, Sechi N, Sfriso A, Sigovini M, Tagliapietra D, Torri- celli P, 2013. The Italian Long-Term Ecosystem Research (LTER-Italy) network: results, opportunities, and chal- lenges for coastal transitional ecosystems. Trans. Wat. Bull. 7:43-63. Quero GM, Perini L, Pesole G, Manzari C, Lionetti C, Bastianini M, Marini M, Luna GM, 2017. Seasonal rather than spatial variability drives planktonic and benthic bacterial diversity in a microtidal lagoon and the adjacent open sea. Mol. Ecol. 26:5961-5973. doi: 10.1111/mec.14363. R Core Team, 2014. A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna. Solidoro C, Bandelj V, Bernardi Aubry F, Camatti E, Ciavatta S, Cossarini G, Facca C, Franzoi F, Libralato S, Melaku Canu D, Pastres R, Pranovi F, Raicevich S, Socal G, Sfriso A, Sigovini M, Tagliapietra D, Torricelli P, 2010. Response of the Venice Lagoon ecosystem to natural and anthro- pogenic pressures over the last 50 years, p. 483-511. In: M. Kennish and H. Paerl (eds.), Coastal lagoons: Critical habi- tats of environmental change. CRC Press, Boca Raton. Schloss PD, Westcott SL, Ryabin T, Hall JR, Hartmann M, Hol- lister EB, Lesniewski RA, Oakley BB, Parks DH, Robinson CJ, Sahl JW, Stres B, Thallinger GG, Van Horn DJ, Weber CF, 2009. Introducing mothur: open-source, platform-inde- pendent, community-supported software for describing and comparing microbial communities. Appl. Environ. Micro- biol. 75:7537-7541. doi: 10.1128/AEM.01541-09 Sfriso A, Facca C, Ceoldo S, Marcomini A, 2005. Recording the occurrence of trophic level changes in the lagoon of Venice over the ‘90s. Environ. Int. 31:993-1001.doi: 10.1016/j.en- vint.2005.05.009 Sfriso A, Buosi A, 2018. Trophic status changes in the Venice Lagoon during the last 40 years. Proceedings XLIX SIBM Congress, Genoa; p. 36-39. Socal G, Bianchi F, Alberighi L, 1999.Effects of thermal pollu- tion and nutrient discharges on a spring phytoplankton bloom in the industrial area of the lagoon of Venice. Vie Mi- lieu 49:19-31. Umgiesser G, Ferrarin C, Cucco A, De Pascalis F, Bellafiore D, Ghezzo M, Bajo M, 2014. Comparative hydrodynamics of 10 Mediterranean lagoons by means of numerical model- ling. J. Geophys. Res. Oceans. 119:2212-2226. doi: 10.1002/2013JC009512. Vadrucci MR, Cabrini M, Basset A, 2007. Biovolume determi- nation of phytoplankton guilds in transitional water ecosys- tems of Mediterranean Ecoregion. Transit. Waters Bull. 2: 83-102. doi: 10.1285/i1825229Xv1n2p83 Wu PF, Li DX, Kong LF, Li YY, Zhang H, Xie ZX, Lin L, Wang DZ, 2020. The diversity and biogeography of microeukary- otes in the euphotic zone of the northwestern Pacific Ocean. Sci. Total Environ. 698:134289 doi: 10.1016/j.scitotenv. 2019.134289 Zingone A, D’Alelio D, Mazzocchi MG, Montresor M, Sarno D, LTER-MC team, 2019. Time series and beyond: multi- faceted plankton research at a marine Mediterranean LTER site. Nat. Cons. 34:273-310. doi: 10.3897/nature conserva- tion.34.30789. No n- co mm er cia l u se on ly