Spatial and Structural Factors Shape Seagrass-Associated Bacterial Communities in Singapore and Peninsular Malaysia

Plant-microbe relationships play critical roles in the functioning and health of terrestrial plants, but little is known about this relationship in marine angiosperms such as seagrasses. Here, we investigated the microbial communities associated with the seagrass Enhalus acoroides throughout Singapore and Peninsular Malaysia. At each sampling location we collected 10 individual and unconnected plants. Each plant was subsequently broken down into leaves, roots, and rhizomes. In addition to living plant parts a sediment sample was taken in close proximity to each. Using high throughput 16S rRNA gene amplicon sequencing we characterised the bacterial communities associated with each plant part and the associated sediment sample. Results indicate geographic structuring of bacterial communities, with a significant pattern of distance decay suggesting dispersal limitation is a contributing factor to the differences we see in bacterial community structure. Bacterial communities can be further differentiated by the function of the collected sample (leaf, root, and rhizome), and we identified a number of microbial indicator species that are associated with each plant part. Further analysis revealed the presence of several microbial taxa that have previously been identified as indicators of “unhealthy” or “stressed” seagrass meadows. This study addresses a current scientific gap related to the characterisation of seagrass microbiomes, and provides a foundation on which future studies can build, particularly those in the Southeast Asian seagrass biodiversity hotspot.

Plant-microbe relationships play critical roles in the functioning and health of terrestrial plants, but little is known about this relationship in marine angiosperms such as seagrasses. Here, we investigated the microbial communities associated with the seagrass Enhalus acoroides throughout Singapore and Peninsular Malaysia. At each sampling location we collected 10 individual and unconnected plants. Each plant was subsequently broken down into leaves, roots, and rhizomes. In addition to living plant parts a sediment sample was taken in close proximity to each. Using high throughput 16S rRNA gene amplicon sequencing we characterised the bacterial communities associated with each plant part and the associated sediment sample. Results indicate geographic structuring of bacterial communities, with a significant pattern of distance decay suggesting dispersal limitation is a contributing factor to the differences we see in bacterial community structure. Bacterial communities can be further differentiated by the function of the collected sample (leaf, root, and rhizome), and we identified a number of microbial indicator species that are associated with each plant part. Further analysis revealed the presence of several microbial taxa that have previously been identified as indicators of "unhealthy" or "stressed" seagrass meadows. This study addresses a current scientific gap related to the characterisation of seagrass microbiomes, and provides a foundation on which future studies can build, particularly those in the Southeast Asian seagrass biodiversity hotspot.

INTRODUCTION
Seagrasses are marine angiosperms whose ancestors made the transition from land back to marine environments between 70 and 100 million years ago (Olsen et al., 2016). Forming large meadows they provide numerous important ecological services, they act as nursery habitats for many marine species, help prevent shoreline erosion through the dissipation of wave energy, reduce human contact with bacterial pathogens and play crucial roles in nutrient cycling and carbon sequestration (Harborne et al., 2006;Waycott et al., 2009;Kannan et al., 2010;Fourqurean et al., 2012;Christianen et al., 2013;Dewsbury et al., 2016;Lamb et al., 2017). It is estimated that seagrasses are responsible for 10% of total global carbon sequestration (Fourqurean et al., 2012). Yet, despite this recognised importance, seagrasses are increasingly threatened by anthropogenic stresses such as coastal development, overfishing, pollution from land based runoff, and climate change (Waycott et al., 2009;Unsworth et al., 2018), and such effects are especially apparent in Southeast Asia (Fortes et al., 2018). Globally, annual seagrass loss is estimated to be approximately 7% per year (Waycott et al., 2009) with declines of approximately 45% reported in localised cases (Yaakub et al., 2014. Plant-microbe associations are critical to the health and functioning of terrestrial plants (Vandenkoornhuyse et al., 2015;Crump et al., 2018), but little is known of the role, diversity and structure of microbial communities in seagrasses and other marine plants (Cúcio et al., 2016;Fahimipour et al., 2017;Crump et al., 2018). However, this is beginning to change as researchers start to examine the microbiomes of marine and coastal plants (Fahimipour et al., 2017;Wainwright et al., 2018Wainwright et al., , 2019bHurtado-McCormick et al., 2019;Lee et al., 2019;Ettinger and Eisen, 2020). Marine plant-associated microbes are involved in sulphide oxidation and reduction, nitrogen fixation, and nutrient cycling. They promote host nutrient uptake and play vital roles in the decomposition of plant matter (Harlin, 1973;Hemminga et al., 1991;Stapel and Hemminga, 1997;Hansen et al., 2000;Crump et al., 2018).
Terrestrial plants show discrete microbial assemblages and communities dependent upon whether leaves, fruits, or roots are examined . Similarly, seagrass associated microbes occupy a number of distinct microenvironments across the plant (e.g., the phyllosphere, and the rhizosphere) (Hurtado-McCormick et al., 2019;Wainwright et al., 2019c). Oxygen levels and the density of inorganic carbon and dissolved organic compounds exuded vary across the seagrass phyllosphere and rhizosphere (Wetzel and Penhale, 1979;Moriarty and Iverson, 1986;Rubio et al., 2017;Ugarelli et al., 2017). Similarly, the roots and rhizomes influence the availability of nutrients in the rhizosphere and create oxygen gradients, and different parts of the seagrass can be submerged, or exposed for differing durations. These differences can promote and lead to the development of bacterial communities specific to each environment or plant part (Ugarelli et al., 2017;Hurtado-McCormick et al., 2019;Wang et al., 2020).
Plants of the same species generally require similar environmental conditions to sustain life, for example, species that are adapted to tropical conditions would not be expected to survive in the tundra, and vice versa. Relatedly, if members of the same species require similar environmental conditions to sustain life, it could be expected that over time plants will select a species-specific consortium of bacteria, and this consortia would likely reflect their requirements for survival and reproduction. Consequently, microbial communities could become homogenised between plant members of the same species, especially when taking into account the high dispersal potential microbes have, or, "everything is everywhere: but the environment selects" (O'Malley, 2008), or in this case the plant selects. Yet, research is frequently showing discrete microbial (fungi and bacteria) communities do exist in different sampling locations within the same species (e.g., seagrasses and mangorves) (Ettinger et al., 2017;Fahimipour et al., 2017;Crump et al., 2018;Wainwright et al., 2019a;Lee et al., 2020).
Southeast Asia is the global epicentre of seagrass biodiversity (Short et al., 2007), while knowledge of seagrass ecosystems is increasing throughout the region, a paucity of research centred on seagrass remains, and considerable knowledge gaps persist . This is especially true of research examining seagrass-associated bacterial communities and microbiomes, which is unfortunate given the increasingly acknowledged vital roles that microbes play in promoting seagrass fitness (Brodersen et al., 2018;Ugarelli et al., 2019Ugarelli et al., , 2017. Understanding these host and microbe relationships is especially pressing in Southeast Asia where seagrass loss is particularly acute, with meadows here suffering from the combined anthropogenic impacts associated with rapid coastal urbanisation and climate change.
With this work we examine the bacterial communities associated with the seagrass Enhalus acoroides throughout Singapore and Peninsular Malaysia and provide foundations on which to build and monitor change. We investigate site-specific associations in an attempt to reveal any patterns of bacterial biogeography, and we profiled above and below ground plant parts to determine whether specific bacterial communities reside in particular structures.

MATERIALS AND METHODS
Ten individual and unconnected Enhalus acoroides plants free of any visible epiphytes were collected at low tide from each location studied (Figure 1 and Supplementary Table 1). To minimise the possible collection of clones, each plant was sampled from a location at least 10 m apart from the last. Collected plants were separated into leaves, roots and rhizomes with a sterile razor blade. Additionally, one sediment sample was collected in close proximity (<1 m) to each plant, with sediment samples taken from approximately 4 cm below the surface. Samples were immediately placed in individual, sealed tubes containing salt saturated CTAB (Hammer et al., 2015), transported on ice to a −20 • C freezer (usually in a guest house). Samples were then transported to the laboratory on ice where they were stored at −80 • C until DNA extraction was performed .
All tissues and sediment samples were disrupted in an Omni Bead Ruptor 24 (Omni International) at 8 m/s for 2 min, and DNA was extracted with a Qiagen DNeasy Powersoil kit following the manufacturer's instructions. PCR amplification targeting the V4 region of the 16S rRNA gene was performed using the bacterial and archaeal primers 515F and 806R (515F -GTG CCA GCM GCC GCG GTA A; 806R -GGA CTA CHV GGG TWT CTA AT). Forward and reverse primers were modified to include Illumina adaptors, a linker and a unique barcode (Caporaso et al., 2011). Each reaction was performed in a total volume of 25 µl, containing 1 µl of undiluted template, 0.1 µl of KAPA 3G Enzyme (Kapa Biosystems, Inc., Wilmington, MA, United States), 0.75 µl of each primer at 10 µM, 12.5 µl KAPA PCR Buffer and water to 25 µl. PCR cycling protocol was 94 • C for 180 s, followed by 35 cycles of 94 • C for 45 s, 50 • C for 60 s and 72 • C for 90 s, with a final extension at 72 • C for 10 min. Negative extraction and PCR controls were included to identify possible contamination issues.
PCR products were visualised on a 1% TBE buffer agarose gel. Normalisation and cleaning of PCR products were performed in SequalPrep normalisation plates (Invitrogen, Frederick, MD, United States) and submitted for sequencing on the Illumina MiSeq platform (600 cycles, V3 chemistry, 300-bp paired end reads) with a 30% PhiX spike (Macrogen Korea).
Sequences were demultiplexed by Macrogen, barcodes and adaptors were removed with Cutadapt (Martin, 2011). Reads were filtered based on quality scores and trimmed using the DADA2 package version 1.14.1 (Callahan et al., 2016) in R version 3.6.2 (R Core Team, 2017). Forward reads were truncated at 250 nucleotides, and reverse reads were truncated at 200 nucleotides. Both forward and reverse reads were filtered to remove any reads with a max EE (expected error) of 2, and reads were additionally truncated at the end of "a good quality sequence" with the parameter truncQ = 2 (see 1 for a detailed explanation of filtering parameters).
The DADA2 algorithm was next used to estimate error rates from all quality-filtered reads and then to merge forward and reverse reads and infer amplicon sequence variants (ASVs). Chimeras were removed with de novo detection. Sequenced extraction negatives were used to identify possible contaminants using the prevalence method implemented in the decontam R package (Davis et al., 2018), and remaining ASVs were assigned taxonomy with the RDP classifier (Cole et al., 2007) 1 https://benjjneb.github.io/dada2/ against a training set based on the Silva v138 16S database (Quast et al., 2013).
Any ASVs assigned to mitochondrial or chloroplast genomes, and those not present in at least 5% of samples were removed. Rarefaction curves were produced using the rarecurve() function implemented in the vegan R package version 2.5-6 (Oksanen et al., 2019). Raw sequence counts were then converted to relative abundance. The Shannon diversity for each sample was calculated, and non-metric multi-dimensional scaling (NMDS) was performed on the Bray-Curtis dissimilarity matrix of samples using the phyloseq R package version 1.30.0 (McMurdie and Holmes, 2013). NMDS plots were generated for all sampled compartments together, and then for each compartment individually (leaf, rhizome, root, sediment). PCoA and network plots were constructed using the phyloseq package. Bar plots of relative abundance were made using ggplot2 version 3.3.1 (Wickham, 2011).
Permutational multivariate ANOVA was performed using the adonis() function of the vegan R package version 2.5-6 (Oksanen et al., 2019). Venn diagrams were created using the VennDiagram R package (Chen and Boutros, 2011). Tests for distance decay of similarity (Mantel test and multiple regression on distance matrices) were performed using the ecodist package. Indicator species were identified for each compartment and each site using the indicspecies R package version 1.7.9 (Cáceres et al., 2012).
All sequences associated with this work have been deposited at the National Center for Biotechnology Information under BioProject ID: PRJNA649070.

RESULTS
After quality control and filtering, 5,539,364 reads were retained for downstream analysis, see Supplementary Table 2 FIGURE 2 | Non-metric multidimensional scaling (NMDS) of bacterial communities based on Bray-Curtis dissimilarity, coloured by location with shapes indicating structure sampled. Plots (A-D), show the leaf, rhizome, root, and sediment samples respectively.
for sequencing statistics of individual samples. Rarefaction curves indicate that sufficient sequencing depth was achieved (Supplementary Figure 1). Bacterial communities appear to be structured by the function of sampled compartments (leaf, rhizome, root, and sediment; Supplementary Figure 2), and these communities can be further differentiated by geographic location (Figure 2). This community structuring is further confirmed by network plots and PCoA plots that show similar patterns of clustering by seagrass part and location (Supplementary Figures 3, 4). Permutational multivariate analysis of variance (PERMANOVA) indicates significant differences in bacterial community according to structure and location (R 2 = 0.183; p = 0.001 and R 2 = 0.092; p = 0.001, respectively; Supplementary Table 3).
Structures that are below ground (i.e., rhizome and root) share more similar bacterial communities than those that are above ground, and biological structures have similar communities and diversity in comparison to non-biological structures (e.g., sediment) (Figure 3 and Supplementary Figure 2). Bacterial communities show a significant pattern of distance decay, with those closer together more similar in community structure (Mantel test: R = 0.128, p = 0.001). This relationship is further supported by multiple regression on distance matrices (MRM) for all samples combined (R 2 = 0.038, p = 0.001) ( Table 1). Bacterial diversity is highest in sediment samples and lowest in the leaves (Figure 3). Samples from Merambong Shoal have the highest median diversity (Shannon), and all locations generally have similar levels of diversity (Supplementary Figure 5). All samples are dominated by members from the phylum Proteobacteria, irrespective of structure and location, and the phylum Fusobacteria has the highest relative abundance in sediment samples from Perhentian (Figure 4). The class Gammaproteobacteria has been found in all compartments, at all locations but tends to be more abundant in living structures (leaf, rhizome, and root samples) (Supplementary Figure 6). Along with Gammaproteobacteria, heatmaps show that Alphaproteobacteria, and Acidimicrobiia and Bacteroidia are common in all compartments. Desulfobacteria and Anaerolineae are more frequently observed in the sediment and below ground structures, and Planctomycetes tends to be most frequently encountered in sediment samples (Figure 5; Supplementary  Figures 6-11).  Amplicon sequence variant (ASV) richness was highest in sediment, followed by root, rhizome, and leaf (1,426, 844, 659, 140, respectively). 64 ASVs are shared between all structures including sediment, and the highest number of shared ASVs (321) are found between rhizome, root, and sediment samples (Figure 6). Indicator species analysis shows leaves are significantly associated with 19 ASVs, rhizomes 16 ASVs, root 48 ASVs and 425 are associated with sediment samples. See Table 2 for the top 5 associations with each structure and Supplementary Table 4 for a complete list. Roots and rhizomes tend to associate with sulphate reducing and oxidising bacteria. The top 5 genera associated with half of the sampling locations have ASVs from the genus Vibrio significantly associated with them (Table 3), see Supplementary Table 5 for full details of all significant associations.

DISCUSSION
Stemming from their high dispersal ability and high abundance, it is a commonly held belief that bacterial communities will show lower levels of spatial variation than larger multicellular organisms (Finlay, 2002;Horner-Devine et al., 2004;Meyer et al., 2018;Bay et al., 2020). Facilitating this expected low spatial variance, microbes are capable of entering extremely lowenergy requiring states where they can persist in dormancy for up to 101.5 Ma. These microbes, when presented with suitable conditions are then able to readily incorporate carbon and nitrogen and go through cell division to rapidly increase in numbers (Bradley et al., 2020;Morono et al., 2020). The ability to enter a dormant state allows bacteria to overcome high environmental heterogeneity, or unfavourable conditions and emerge when conditions are suitable (Jones and Lennon, 2010). This idea suggests that while microbes could be everywhere, they might not actually be playing biological or ecological roles in the ecosystem at the particular point in time sampling took place. Molecular methods would still detect their presence irrespective of whether they are active or dormant, likewise, relic DNA (that from dead microorganisms) can persist for years (Cangelosi and Meschke, 2014;Carini et al., 2017). The presence of dormant microbes and relic DNA can obscure biologically relevant patterns and phenomena and suggest an absence of geographically localised community structure (Finlay, 2002;Carini et al., 2017).
However, work is now challenging the notion that variance is low in microbial communities, and congruent with other research examining marine microbiomes (Cúcio et al., 2016;Crump et al., 2018;Bay et al., 2020;Osman et al., 2020;Tan et al., 2020;Wainwright et al., 2020) we show that bacterial communities can be significantly different between sampling location and structure examined, and this spatial variance is particularly strong in sediment samples in comparison to living seagrass structures. This is consistent with previously proposed hypotheses suggesting that habitats offered by living plant organs (i.e., leaf, fruit, etc.) act as a biological filter allowing the plant to exert some degree of control over the constituents of their microbial communities (Goldmann et al., 2016;Jones et al., 2019;Wang et al., 2020). Consequently, the microbial communities associated with living compartments should be more similar to each other over wider geographic scales, whereas soil and sediment associated communities are not subjected to any filtering, therefore we see higher community variance between sample sites in these non-living sample types, and this could allow the development of more unique communities. This idea is supported by the results of our ordinations, MRM and Mantel tests. Ordinations show tight and distinct clustering in sediment samples, and a greater degree of community overlap in living structures, but there are still significant differences in community structure between living parts. MRM and Mantel tests confirm this structuring and the strongest patterns of distance decay of similarity in bacterial communities is found in sediment samples. Further supporting this idea of filtering we see the highest microbial diversity in sediment samples and much lower, but similar level of diversity in all living samples. This likely reflects the specific requirements that each living structure requires for correct functioning, and similar to our work, numerous other studies report highly diverse microbial communities in soil and sediment samples (Serna-Chavez et al., 2013;Fierer, 2017).
Our analysis of indicator species shows a number of microbes significantly associated with each structure and location, for FIGURE 4 | Stacked bar plots of relative bacterial abundance. All samples of a specific type, from one location combined. example, rhizome and root tissues both associate with bacteria that are involved in sulphide-oxidation. Sulphur containing compounds are highly toxic to seagrasses and cause a reduction in photosynthetic performance and block aerobic respiration by inhibition of the mitochondrial cytochrome c oxidase enzyme (Lee et al., 1999). Under anoxic conditions sulphide compounds (H 2 S, HS − , and S2 − ), especially hydrogen sulphide, a common metabolic poison can accumulate to levels that are toxic (Goodman et al., 1995). These sulphide containing compounds enter root tissue where they inhibit enzyme performance and functioning (Lee et al., 1999), and field experiments show that sulphide exposure results in stunted growth (King et al., 1982). The presence, and significant associations of microbes that can reduce the availability of toxic compounds in these structures is a likely adaption to the anoxic conditions that exist in water logged sediment environments. Leaves were significantly associated with bacteria from the genus Marinomonas. This genus has previously been found associated with the leaves of the seagrass Posidonia oceanica (Garcias-Bonet et al., 2012;Tarquinio et al., 2019), and plays a role in the production of oxidative enzymes such as melanin (Sanchez-Amat et al., 2010). Melanin is a free radical scavenger that helps prevent cellular damage that can be caused by oxidative stress at high light intensities and the resultant excess excitation energy (Mittler, 2002;Costa et al., 2015;Tarquinio et al., 2019). Microbes that minimise oxidative stress are particularly useful in the leaves of seagrasses from tropical regions where exposure to direct and intense sunlight at low tide is frequent. Taxa from the genus Marinomonas are also recognised as Plant Growth Promoting Bacteria (PGPB), which assist and increase growth and development in macroalgae (Singh et al., 2011), and it is suggested that the same PGBPs have profound effects and enhance the development rate of seagrass leaves (Celdrán et al., 2012). Similarly, leaves were significantly associated with the genera Labrenzia which also helps promote plant growth and has antimicrobial properties (Martin et al., 2020).
Microbial indicators have been proposed as sensitive, and effective markers of environmental perturbation and stress, and the development/identification of these markers and their associations can aid in ecosystem management and restoration programmes (Glasl et al., 2019;Wainwright et al., 2019a;Martin et al., 2020). Previous studies have identified a number of   Thiogranum Rubripirellula bacterial genera that are indicative of "healthy" and "stressed" seagrass states. For example, methylotrophic bacteria, iron cycling bacteria, and N 2 fixing bacteria are generally more abundant and associated with seagrass meadows identified as "healthy." Conversely, sulphur-cycling bacteria (sulphideoxidising and sulphate-reducing) are more frequently observed in "stressed" meadows (Martin et al., 2020), and the genus Vibrio has been implicated as a causative agent of disease in many marine species (Colwell and Grimes, 1984;Jayasree et al., 2006;Vezzulli et al., 2016;Weynberg et al., 2016;Rubio-Portillo et al., 2018). Several studies have putatively identified Vibrio spp. as bacterial pathogens in seagrass meadows, and Vibrio spp. tend to show increased abundance in more impacted sites (Liu et al., 2018;Tarquinio et al., 2019). Indicative of the anthropogenic stresses associated with coastal urbanisation and increasing coastal development throughout the region, several of our sample sites have significant associations with Vibrio spp. and sulphur cycling bacteria. We did not detect any of the bacterial genera that have been identified as indicators of "healthy" seagrass meadows, and given the conservation challenges that seagrasses in the region face, this is not entirely surprising. However, caution should be exercised when interpreting what bacteria are indicative of "healthy" seagrass meadows, especially when comparing study results from different countries, particularly given the complexities and dynamics of microbial systems (Archer et al., 2019). This study is the first of this nature from the region, and the composition of the microbial constituents in a healthy seagrass meadow here could be very different to those in other regions. Nevertheless, a number of the bacteria we identified are associated with degraded seagrass meadows and it would be wise to consider these findings in seagrass conservation schemes. Doing so could facilitate the taking of Desulfofustis steps to mitigate known stressors such as increased sedimentation resulting from deforestation, particularly as these stresses have already been documented to increase the prevalence of potential seagrass pathogens (Liu et al., 2018). As has previously been suggested (Santos et al., 2011;O'Callaghan, 2016;Zahn and Amend, 2017;Lee et al., 2019;Tarquinio et al., 2019;Vanwonterghem and Webster, 2020), understanding microbial community dynamics and structure is likely a key determinant of restoration and conservation success. Unfortunately, the majority of seagrass restoration projects end in failure (Bayraktarov et al., 2016;van Katwijk et al., 2016) and a variety of reasons are thought to contribute to these unfavourable outcomes, amongst them -poor site location, techniques used, or ongoing human stressors (Fonseca, 2011). To our knowledge, and despite the widespread manipulation of microbiomes in terrestrial systems to achieve positive outcomes (O'Callaghan, 2016;Zahn and Amend, 2017), this has not been attempted in seagrass restoration projects. Considering the microbial constituents in seagrass restoration is likely even more important when using transplants from other regions, or transplants grown in ex situ nurseries, especially, as we show that microbial communities can differ over comparatively small spatial scales (<10 km). If microbes are not considered it is possible that transplants could be maladapted to their new environments, or suffer from transplantation shock (Wang et al., 2021) which could be a contributing factors leading to restoration failure (Vogel et al., 2021).
With this work we characterise the bacterial communities associated with the seagrass Enhalus acoroides throughout Singapore and Peninsular Malaysia, providing a framework on which future studies can build and monitor bacterial community change and composition.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm. nih.gov/, PRJNA649070.

AUTHOR CONTRIBUTIONS
All authors listed have made a substantial, direct and intellectual contribution to the work, and approved it for publication.

FUNDING
This study was funded by the National Research Foundation, Prime Minister's Office, Singapore under its Marine Science R&D Programme (MSRDP-P03) and the Mandai Nature Fund. Student support was provided by the Summer Research Programme, administered by the Yale-NUS College Centre for International & Professional Experience (CIPE) and co-funded by the Yale-NUS Dean of Faculty Office. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

ACKNOWLEDGMENTS
All samples were collected under permit number NP/RP16-156 issued by the National Parks Board of Singapore. Collections from Malaysia were made under permit JTLM 630-7Jld.9(9).