Eelgrass Leaf Surface Microbiomes Are Locally Variable and Highly Correlated with Epibiotic Eukaryotes

Eelgrass (Zostera marina) is a marine foundation species essential for coastal ecosystem services around the northern hemisphere. Like all macroscopic organisms, it possesses a microbiome (here defined as an associated prokaryotic community) which may play critical roles in modulating the interaction of eelgrass with its environment. For example, its leaf surface microbiome could inhibit or attract eukaryotic epibionts which may overgrow the eelgrass leading to reduced primary productivity and subsequent eelgrass meadow decline. We used amplicon sequencing of the 16S and 18S rRNA genes of prokaryotes and eukaryotes to assess the leaf surface microbiome (prokaryotes) as well as eukaryotic epibionts in- and outside lagoons on the German Baltic Sea coast. Prokaryote microbiomes varied substantially both between sites inside lagoons and between open coastal and lagoon sites. Water depth, leaf area and biofilm chlorophyll a concentration explained a large amount of variation in both prokaryotic and eukaryotic community composition. The prokaryotic microbiome and eukaryotic epibiont communities were highly correlated, and network analysis revealed disproportionate co-occurrence between a limited number of eukaryotic taxa and several bacterial taxa. This suggests that eelgrass leaf surfaces are home to a mosaic of microbiomes of several epibiotic eukaryotes, in addition to the microbiome of the eelgrass itself. Our findings thereby underline that eukaryotic diversity should be taken into account in order to explain prokaryotic microbiome assembly and dynamics in aquatic environments.

Eelgrass (Zostera marina) is a marine foundation species essential for coastal ecosystem services around the northern hemisphere. Like all macroscopic organisms, it possesses a microbiome (here defined as an associated prokaryotic community) which may play critical roles in modulating the interaction of eelgrass with its environment. For example, its leaf surface microbiome could inhibit or attract eukaryotic epibionts which may overgrow the eelgrass leading to reduced primary productivity and subsequent eelgrass meadow decline. We used amplicon sequencing of the 16S and 18S rRNA genes of prokaryotes and eukaryotes to assess the leaf surface microbiome (prokaryotes) as well as eukaryotic epibionts in-and outside lagoons on the German Baltic Sea coast. Prokaryote microbiomes varied substantially both between sites inside lagoons and between open coastal and lagoon sites. Water depth, leaf area and biofilm chlorophyll a concentration explained a large amount of variation in both prokaryotic and eukaryotic community composition. The prokaryotic microbiome and eukaryotic epibiont communities were highly correlated, and network analysis revealed disproportionate co-occurrence between a limited number of eukaryotic taxa and several bacterial taxa. This suggests that eelgrass leaf surfaces are home to a mosaic of microbiomes of several epibiotic eukaryotes, in addition to the microbiome of the eelgrass itself. Our findings thereby underline that eukaryotic diversity should be taken into account in order to explain prokaryotic microbiome assembly and dynamics in aquatic environments.

INTRODUCTION
Seagrasses are aquatic flowering plants that form underwater meadows critically important for coastal ecosystems around the world. Seagrass meadows are nursing grounds for juvenile fish which can hide and forage between the seagrass leaves, sediments are stabilized by the seagrass roots and the biomass of seagrass and associated organisms sequester carbon with implications for climate change mitigation (Fourqurean et al., 2012). The seagrass eelgrass (Zostera marina) is an important foundation species along soft-bottom coasts in the northern hemisphere. Eelgrass meadows are biodiversity hotspots, providing a home to a myriad of invertebrates, algae, fish and microorganisms. Due to its ability to tolerate low and fluctuating salinity levels, eelgrass is common in estuaries as well as in the biggest brackish water habitat in the world, the Baltic Sea. However, significant declines in the depth limit and areal cover of eelgrass meadows in the Baltic Sea have been observed since several decades (Boström et al., 2014). This echoes the distressing situation for seagrass ecosystems around the world, which experience threats by human activities and global change (Orth et al., 2006).
The mechanisms that are responsible for eelgrass meadow decline appear to depend on many different factors. Historically, outbreaks of the pathogenic protist Labyrinthula zosterae ("wasting disease") has caused catastrophic die-off of eelgrass (Muehlstein et al., 1991), although its importance as a pathogen under current conditions in Europe seems limited (Brakel et al., 2014). Instead, factors such as water clarity, eutrophication, grazing pressure, and interspecific competition have been identified as culprits for eelgrass productivity decline (Baden et al., 2010;Duffy et al., 2015). Especially important is the competition between eelgrass and algae that grow within meadows, often as epiphytes on the eelgrass itself. In some cases, eelgrass becomes extensively covered with a mixture of algae, bacteria and sessile animals, which shade it and inhibit transport of solutes (Brodersen et al., 2015) and may over time cause eelgrass meadows to degrade. Recent research suggests that relative success of eelgrass and associated algae is determined by complex interactions between biotic processes such as grazing on both algae and eelgrass by invertebrates, and environmental factors such as nutrient concentrations and temperature (Eklöf et al., 2012;Alsterberg et al., 2013). However, these mechanisms are not well understood and their complexity necessitates a holistic ecosystem approach taking into account several organism groups that inhabit seagrass meadows and their interactions to be resolved (Boström et al., 2014;Maxwell et al., 2016).
An overlooked group of organisms in eelgrass meadows is the prokaryotes (bacteria and archaea), which cover eelgrass leaves and roots forming biofilms, are responsible for degradation of eelgrass detritus and are also associated to all other organisms in eelgrass beds including algae and other epiphytes. Prokaryotes are the first colonizers on new seagrass leaves and thus initiate a successional process that may end with severe epiphytic overgrowth. They are also an important food source for microbial eukaryotic and invertebrate grazers and may therefore support a substantial part of the foodweb found within eelgrass meadows. Research on eelgrass leaf prokaryotic communities is very limited, but early work has determined that prokaryotic abundance and productivity on eelgrass leaves vary during the year, with a peak in early autumn (Törnblom and Søndergaard, 1999), and that eelgrass leaves have an associated prokaryotic community (the "microbiome"), with only some compositional overlap with other aquatic macrophytes (Crump and Koch, 2008).
There is currently little known about what factors influence eelgrass microbiomes, or what functions the prokaryotes have on eelgrass leaves. Bacteria and archaea are likely to play a fundamental role in the competition of eelgrass and epibionts such as algae. For example, certain bacteria may inhibit the attachment of algal spores while others may promote further colonization (Celdrán et al., 2012;Mieszkin et al., 2013). Conversely, epibiotic eukaryotes including algae may in turn shape the prokaryotic communities on eelgrass leaves by release of chemical attractants or deterrents (Steinberg and de Nys, 2002) or through selective grazing (Huws et al., 2005), for example. Biotic and abiotic environmental factors as well as host-related factors such as eelgrass productivity and genotype are also likely to shape the eelgrass microbiome, and thereby modulate its function.
In this study, we aimed to obtain a first view onto eelgrass leaf microbiomes by investigating the community composition of prokaryotic and eukaryotic epibiotic organisms in relation to abiotic and biotic environmental variables. We sampled eelgrass in semi-sheltered lagoons and along exposed open shorelines around the Island of Hiddensee on the eastern German Baltic Sea coast and used high-throughput Illumina amplicon sequencing of the 16S and 18S rRNA genes. We hypothesized that (1) both prokaryotic and eukaryotic communities would vary substantially between lagoon and open coast sites due to different abiotic conditions and that (2) prokaryotic microbiome composition would be influenced by the composition of epibiotic eukaryotes.

Study Area
We chose to perform our survey in the area around the island of Hiddensee on the German Baltic Sea coast as it offers variable yet locally representative environments colonized by eelgrass meadows. The area is characterized by semienclosed, shallow lagoons (German: "Bodden") which differ from the open coast with respect to depth distribution, salinity and nutrient concentrations. The surveyed lagoons (the Vitter Bodden and Schaproder Bodden) feature somewhat lower average salinity (8.8 ± 1.0 , range 6.5-13.4 ) compared to the open coast (Libben bay, 9.4 ± 1.6 , range 6.8-15.9 ). Nutrient levels are elevated in lagoons (total N: 38.2 ± 12.3, total P: 1.2 ± 0.62 µmol l −1 ) compared to open coast waters (total N 19.9 ± 4.3, total P: 0.91 ± 0.4 µmol l −1 ) mainly due to agricultural runoff and other human activities (Schiewer, 2008). Salinity and nutrient values are yearly averages from monitoring data 2005 -2014 (Landesamt für Umwelt, Naturschutz und Geologie Mecklenburg-Vorpommern, unpublished data).

Sampling
Zostera marina shoots were sampled via scuba diving at seven sites around the island of Hiddensee (Figure 1). At every site, five replicate shoots (above-ground parts) were collected within an area of approximately 1 m 2 . Shoots were kept cool and in clean plastic bags until sample processing which was carried out within a few hours. Environmental variables such as water depth and parameters relating to the surrounding macrophyte vegetation (eelgrass and co-occurring macroalgae and aquatic plants) were determined as part of a simultaneous macrophyte inventory (Bühler, 2016).

Leaf Biofilm Sample Processing
The sampled eelgrass shoots were rinsed with sterile filtered seawater and handled under sterile conditions. Leaf areas with heavy macroscopic epibiosis and visibly degrading tissue were removed before biofilm sampling. Epibiotic biofilm from the entire leaf surface of each shoot was scraped off using sterile cotton swabs. The scraped off material was suspended in sterilefiltered seawater (10 ml) and kept cool. The suspensions were aliquoted into separate tubes and centrifuged (10000 × g, 10 min, 4 • C) to pellet biofilm material. Pellets were frozen at −20 • C until DNA and chlorophyll a extraction. Chlorophyll a was extracted overnight in acetone at +4 • C and absorbance at 665 and 750 nm (A 665 and A 750 ) was measured spectrophotometrically (Genesys 20, Thermo Fisher Scientific). Biofilm chlorophyll a content (C chla ) was calculated as C chla = 11.4 * (A 665 -A 750 ) and was expressed in µg chlorophyll a cm −2 of biofilm.

Sequence Processing
Sequences clipped from adaptor and primer sequence remains were processed using the DADA2 package in R (version 1.2.0) (Callahan et al., 2016;R Development Core Team, 2017). Briefly, forward and reverse Illumina reads were truncated to 200 bp, filtered (maxEE = 2, truncQ = 2), dereplicated and error rates were estimated using the maximum possible error estimate from the data as an initial guess. Sample sequences were inferred, and paired forward and reverse reads were merged. Chimeric sequences were removed using the removeBimeraDenovo function. The resulting unique sequence variants (analogous to operational taxonomic units) were used to construct a table containing relative abundances of sequence variants across all samples. Sequence variants were taxonomically classified using a lowest common ancestor approach (Lanzén et al., 2012) based on the Silva database (Pruesse et al., 2007).

Statistical Analyses
Multivariate statistical analyses were carried out using functions in the vegan package [version 2.4-1, (Oksanen et al., 2016) To visualize similarities in sequence variant composition (community composition) between sampling sites, non-metric multidimensional scaling (vegan function metaMDS) was performed on Hellinger-transformed sequence variant counts using Bray-Curtis distance. To explain the variation in community composition in response to the environment, abiotic and biotic variables were selected based on their explanatory power in PERMANOVA tests (vegan function adonis, variables were transformed (log e ) to achieve normaldistribution if needed). A final model including the two most relevant environmental variables was selected and distance-based redundancy analysis ordination (vegan function capscale, Bray-Curtis distance) constrained to these environmental variables was used to visualize variation.
Sequence variant richness was calculated by rarefying read counts to the lowest number of reads in a sample of each dataset (16S: 5861, 18S: 25780). Evenness was calculated as J = H/logS, where H is the Shannon diversity index and S is rarefied richness (Pielou, 1977). To identify sequence variants that were differently abundant in lagoon sites and open coast sites, differential abundance analysis as implemented in the DEseq2 R package (Love et al., 2014).
A partial mantel test was used to assess correlation between prokaryotic and eukaryotic community composition while correcting for environmental variables [vegan function mantel.partial, Bray-Curtis distance for community data, euclidean distance for environmental variables, (Legendre and Legendre, 2012)]. To analyze pairwise correlations between 16S and 18S sequence variants, co-occurrence networks were calculated in R using Pearson correlation. P-values were adjusted for multiple testing using the Benjamini-Hochberg method (Benjamini and Hochberg, 1995). Only correlations with adjusted p-values < 0.01 were treated as significant and included in the network. The network was plotted and edited for clarity (symbols, colors etc.) in Cytoscape version 3.3.0 (Shannon et al., 2003). (Bühler, 2016), all open coast sites where eelgrass was found were deeper than lagoon sites. Eelgrass biomass (leaf surface area and dry weight) varied widely among sampling sites ( Table 1). Leaf surface biofilms were visually different in color and density (results not shown) and biofilm chlorophyll a content varied substantially among sampling sites ( Table 1).

Due to the different depth distributions of eelgrass on open coasts and in lagoons
Illumina amplicon sequencing of 16S and 18S rRNA fragments resulted in 1.0 million 16S reads and 1.9 million 18S reads with an average of 29564 16S reads (min = 5861, max = 67456) and an average of 56886 18S reads per sample (min = 25780 max = 126938). One sample from site SBL was excluded from all analyses because of low 16S read count (861 reads). In total, 4409 prokaryotic and 1946 eukaryotic sequence variants (analogous to operational taxonomic units) were identified. Among prokaryotes, Bacteria were overwhelmingly abundant over Archaea (99.2% of sequence variants and 99.9% of reads). Due to this, 16S sequence variants are referred to as "bacteria" from here on.
Bacterial community composition varied substantially among sites (PERMANOVA R 2 = 0.82, p < 0.001), and the main variation was between lagoon sites and open coast sites (PERMANOVA R 2 = 0.42, p < 0.001, Figure 2A). Eukaryotic community composition also varied among individual sites (PERMANOVA R 2 = 0.83, p < 0.001) and varied strongly between lagoon and open coast sites (PERMANOVA R 2 = 0.48, p < 0.001, Figure 2B). Water depth, eelgrass shoot leaf area and biofilm chlorophyll a content were the strongest environmental predictors of community variation in both bacteria and eukaryotes, with depth explaining about 30%, leaf area 11% and chlorophyll a 10% of total variation (PERMANOVA p < 0.001). There was also a significant interaction term between depth and chlorophyll a content for both bacterial and eukaryote communities (PERMANOVA, p < 0.01), indicating that chlorophyll a affects communities differently depending on depth, or vice-versa. The variation explained by depth and chlorophyll a content is illustrated in Figures 2C,D. Richness of bacterial sequence variants was positively correlated with eelgrass dry weight at the sampling site (R 2 = 0.69, p < 0.001, Figure 3A) as well as with eelgrass leaf area, although somewhat more weakly (R 2 = 0.29, p < 0.001). Eukaryotic sequence variant richness did not correlate with eelgrass dry weight or other environmental gradients, yet was significantly higher in lagoon sites compared to open coast sites (ANOVA F = 208.6, p < 0.001, Figure 3B). Richness and evenness were correlated for both bacteria (R 2 = 0.34, p < 0.001) and eukaryotes (R 2 = 0.84, p < 0.001). A differential abundance analysis revealed sequence variants that differed significantly (adjusted p < 0.01) in relative abundance between lagoon and outer coast sites (Figure 4). On a broad taxonomical level, several bacterial sequence variants from Cyanobacteria, Planctomycetes, Gemmatimonadetes, and Acidobacteria were overrepresented in lagoon sites ( Figure 4A). Among eukaryotic taxa, one highly abundant individual sequence variant, classified as Hydrozoa (Metazoa: Cnidaria) was significantly overrepresented in open coast sites, while several other metazoan (e.g., classified as Dorvillea and Rotifera) and diatom sequence variants (Cocconeidaceae) were overrepresented in lagoon sites. Fourteen bacterial and 3 eukaryotic sequence variants were detected in all samples, and were defined as the "core community" of eelgrass leaf biofilms for the purpose of this study ( Table 2).
A partial Mantel test, correcting for all measured environmental variables, confirmed a strong correlation Leaf area and chlorophyll a values represent the mean of five replicates ± 1 standard deviation. Water depth is the depth at which eelgrass was collected.  between bacterial and eukaryote community composition (partial mantel r = 0.92, p = 0.001), which was suggested by the similar clustering in the nMDS ordinations Figures 2A,B.
To further explore correlations between individual bacterial and eukaryotic sequence variants, a co-occurrence network was constructed (Figure 5). The network features several prominent sequence variants that correlate with many other taxa and in some cases form distinct network modules. In particular, a sequence variant classified as Hydrozoa (Metazoa: Cnidaria) displayed both positive (co-occurrence) and negative (mutual exclusion) correlations with several bacterial taxa. Other prominent eukaryotic taxa were classified as Mytilus sp. (mussel), Dorvillea sp. (polychaete) and Bacillariophyta (diatoms). In some cases, network modules centered around bacteria, such as Bacteroidetes, Proteobacteria, and Cyanobacteria.

Eelgrass Bacterial Microbiomes Are Locally Variable
This study represents the first effort to characterize both bacterial and eukaryotic epibiotic communities on eelgrass leaves and to relate their communities to potentially important environmental drivers and to each other. We found that the local abiotic environment, in the form of water depth, explained a great deal of variation for both communities. Water depth was one of the main differences between lagoon and open coast sites which likely in part explains the clear separation observed between these sites. However, water depth encompasses separate physical factors such as light penetration and wave shear forces which are difficult to further separate in this study. In addition, several other factors such as mean salinity, nutrient levels and water temperature differ on seasonal time scales in lagoons compared to the open coast (Schiewer, 2008;Blindow et al., 2016) and could also contribute to this variation. Of the biotic environmental variables measured, eelgrass shoot leaf area (a measure of the size of each shoot) and biofilm chlorophyll a content explained the most variation in community composition of both bacteria and eukaryotes. The explanatory power of eelgrass leaf area, which is related to eelgrass productivity (Echavarría-Heras et al., 2010), suggests an influence of host physiology on epibiotic biofilm composition. The amount of variation explained by chlorophyll a content, a proxy for the contribution of photosynthetic organisms (algae and cyanobacteria) in the biofilms, further indicates that interactions within the biofilms, for example between eukaryotic algae and bacteria play an important role in shaping community composition.
Richness of bacterial sequence variants was correlated with leaf area and with eelgrass dry weight at the sampling sites. Larger leaf area corresponds directly to a larger area of biofilm sampled, since biofilm from entire eelgrass shoots was removed. This observed richness-area relationship agrees with similar biodiversity-area relationships in other microbial habitats (reviewed in Martiny et al., 2006) and among macroscopic organisms. However, an alternative explanation for this pattern relates to the classical notion that energy supply limits biodiversity (Currie, 1991). Since leaf area is a proxy for eelgrass productivity (Echavarría-Heras et al., 2010), its correlation with bacterial diversity could indicate that higher resource supply in the form of for example organic exudates promotes leaf surface bacterial diversity. In contrast, eukaryotic biodiversity (measured as sequence variant richness) seemed to be controlled by factors specific to lagoon vs. open coast habitats, rather than area.
Despite the small geographic area surveyed in this study, variability in community composition of both bacteria and eukaryotes was substantial among sites, which shows that eelgrass leaf biofilms are dynamic habitats on a regional scale. This observation is consistent with, and complementary to, recent results that revealed that eelgrass leaf surface microbiomes were more variable than root-associated microbiomes on a global scale (Fahimipour et al., 2017), as well as within an individual eelgrass patch (Ettinger et al., 2017). Nevertheless, we detected several sequence variants in all samples analyzed, supporting the possible existence of a core microbiome of eelgrass leaves. Several core sequence variants were classified to the family Rhodobacteraceae (Alphaproteobacteria), which are often surfaceattached, chemoorganotrophic bacteria common in the marine environment (Dang and Lovell, 2002). Interestingly, this family was also consistently detected in the rhizospheres of eelgrass and its close relative Z. noltii (Cúcio et al., 2016;Ettinger et al., 2017) as well as on leaves and below-ground parts of the seagrass Halophila stipulacea in the red sea (Mejia et al., 2016) indicating that their members indeed play important roles in association with various seagrasses. The observation of a Methylophilaceae (Betaproteobacteria) sequence variant as an abundant core community member agrees with detection of this taxon on eelgrass and other co-occurring macrophytes in the brackish-water Chesapeake Bay (Crump and Koch, 2008) and suggests importance of one-carbon metabolism on eelgrass leaf surfaces (Lapidus et al., 2011).

Eelgrass Bacterial Microbiomes Correlate with Eukaryotic Epibiont Communities
The strong correlation between the bacterial and eukaryotic communities observed in this study has at least two possible causes. First, the abiotic and biotic environment experienced by both communities could influence them in similar ways. This Ectocarpales explanation is partially supported by the similar amounts of variation explained by water depth, leaf area and chlorophyll a content for both communities. Second, bacterial and eukaryotic taxa may interact directly, resulting in a bacterial community that is driven by the composition of the eukaryotic community, or vice versa. The topology of the co-occurrence network (Figure 5) suggests some importance of the latter scenario, where certain eukaryotic taxa (e.g., classified as Hydrozoa, Mytilus sp. and Cocconeidaceae) feature a disproportionate number of significant correlations with bacterial taxa. This indicates that bacterial communities in mature biofilms on eelgrass leaves are to some extent a reflection of the combined microbiomes of several epibiotic eukaryotes such as hydrozoans, mussels and diatoms in addition to any eelgrass-specific microbiome that may exist. Conversely, some bacterial taxa correlate with several eukaryotic taxa. For example, a cyanobacterial sequence variant classified as Rivularia sp. correlates with several ciliate taxa which may for example be attracted to this cyanobacterium as grazers.
In a previous study, bacterial and eukaryotic communities of the tropical seagrass Enhalus acoroides were found to respond similarly to environmental CO 2 exposure, suggesting that seagrass bacterial and eukaryotic epibiotic communities are highly correlated also in other seagrass species and environments. However, potential direct interactions between bacterial and eukaryotic taxa were not addressed in this study (Hassenrück et al., 2015).

Outlook
To date, many studies addressing the microbiomes of marine organisms have focused exclusively on prokaryote communities (see e.g., Hentschel et al., 2012;Wahl et al., 2012;Bourne et al., 2016 and references therein). In contrast, free-living microbial eukaryotes are receiving increasing attention and molecular community analysis has led to several recent discoveries that emphasize their underexplored diversity and crucial roles in the environment (Liu et al., 2009;Lima-Mendez et al., 2015). It is likely that insights of a similar caliber can be gained by applying molecular tools to eukaryotes associated to macroscopic hosts (Andersen et al., 2013), beyond relatively well-described eukaryotic symbionts such as zooxanthella in corals (Muller-Parker et al., 2015). We have shown that the bacterial microbiome is strongly correlated with epibiotic eukaryote communities, including microbial and macroscopic taxa, indicating that interactions between eukaryotes and bacteria drive community assembly on eelgrass surfaces to some extent. By ignoring eukaryotes in bacterial microbiome studies, a large part of the variation in bacterial communities could remain unexplained or misinterpreted. Indeed, the high diversity of both eukaryotic primary producers such as diatoms and heterotrophic taxa such as ciliates, rhizarians and many macroscopic and microscopic animals paint a picture of complex microbial ecosystems encompassing several trophic levels on eelgrass leaves. In these "microbial jungles, " bacteria are bound to play important roles which may be especially relevant in interaction with their eukaryote neighbors. In addition to considering the understudied eukaryotes in microbiome studies, gaining a mechanistic understanding FIGURE 5 | Co-occurrence of bacterial and eukaryotic rRNA gene sequence variants presented as a network. Solid lines (network edges) indicate significant positive Pearson correlations (r > 0.8, p < 0.01) between bacterial (circles) and eukaryotic (triangles) rRNA gene sequence variants (i.e., taxa, network nodes) while dashed lines indicate significant negative correlations (r > -0.8, p < 0.01). Nodes are sized according to the number of correlating neighbors (degree) and color-coded according to taxonomic affiliation. Selected prominent nodes (same as in Figure 4) are further labeled with more detailed taxonomic classifications.
into microbiome assembly requires means to infer causative relationships between microbiome dynamics and abiotic and biotic factors (Widder et al., 2016). The present study included a limited number of samples from a small geographic area, providing a valuable first assessment of environmental drivers of community assembly and potential microbial interactions. Further insight could be gained by manipulative experiments, whereby for example, seagrass surfaces are exposed to colonization by microorganisms under controlled conditions. In addition, large observational microbiome datasets from a range of eelgrass meadows encompassing wider environmental gradients would allow more statistical power to tease apart different drivers of community assembly. Modeling approaches such as structural equation modeling (Pugesek et al., 2003) could allow assessment of causative relationships between the environment, the seagrass host and its microbiome. A major challenge lies in linking microbiome dynamics to processes that are relevant on the ecosystem level.
Manipulative experiments that operate on an ecosystem scale, i.e., mesocosm experiments mimicking seagrass meadows, offer a promising approach to manipulate the environment while tracking microbiome dynamics and their links to critical ecosystem fluxes such as primary production and respiration. Ultimately, understanding assembly and function of seagrass microbiomes could contribute to tackling some of the challenges faced by seagrass meadow ecosystems under global climate change and other anthropogenic stressors.

AUTHOR CONTRIBUTIONS
MB conceived the study, participated in the field work, processed the eelgrass samples, led the laboratory work, analyzed the data and wrote the paper. ABü planned and executed the scuba diving, gathered the macrophyte data and prepared Figure 1. ABr performed molecular laboratory work and contributed to data analysis. SD participated in planning and execution of the field work and drove the research vessel. HS planned the fieldwork. IB planned and supervised the fieldwork and macrophyte data collection. All authors edited the manuscript.

FUNDING
The field sampling campaign was funded through a grant from the German Federal Ministry of Education and Research (Project BACOSA, 03F0665C) awarded to HS and IB. Microbial community analysis was funded through the University of Greifswald.