High-throughput sequencing characterizes intertidal meiofaunal communities in Northern Gulf of Mexico (Dauphin Island and Mobile Bay, Alabama).
Meiofauna are generally defined as animals 45 [micro]m to 1 mm that live between sediment grains (Mare, 1942; Giere, 2009). These organisms are a key component of food webs, nutrient exchange, and sediment stability (Snelgrove, 1997). More specifically, meiofauna serve as a trophic linkage between bacteria and macrofauna (Coull, 1999) as well as play an important role in transport between the sediment and water column (Aller and Aller, 1992; Coull, 1999). Despite their ecosystem importance, meiofauna are often overlooked and not easily studied because sieving and using microscopy to taxonomically identify individuals is extremely time consuming and usually requires experts to properly identify organisms. Furthermore, because expertise is needed, researchers typically focus on a few taxonomic groups within meiofauna, limiting the ability of studies to be extrapolated across the community. High-throughput sequencing approaches along with bioinformatic tools provide the ability to examine the composition of communities as a whole more quickly than traditional methods, allowing for a better understanding of community composition and variation. Even with the usefulness of this technology, only a handful of studies have implemented high-throughput genomic approaches to assess meiofaunal diversity (Fonseca et al, 2010; Creer et al., 2010; Bik et al., 2012a, b, c).
In recent years, interest in the Gulf of Mexico (GOM) ecosystem has grown greatly in light of the Deepwater Horizon Oil Spill (DWHOS). Unfortunately, current knowledge of meiofauna community composition and variation within the GOM is limited, especially in intertidal systems. The few studies that have examined subtidal meiofauna community composition in GOM focused either on deep sea (Pequegnat et al., 1990; Baguley et al., 2006, 2008; Montagna et al., 2013), on continental shelf locations (Montagna and Harper, 1996; Escobar-Briones and Soto, 1997; Escobar et al., 1997; Landers et al., 2012), or within Texas estuaries (Montagna and Kalke, 1992; Montagna et al., 2002). Nematodes and copepods were reported as the two dominant taxa within communities (Montagna et al., 2002; Baguley et al., 2006, 2008; Landers et al., 2012). Additionally, variation in meiofauna abundance has been related to water depth and longitudinal locations (Pequegnat et al., 1990; Baguley et al., 2006; Landers et al., 2012). These studies provided baseline knowledge of community composition and allowed for exploration of the impact of the DWHOS, which greatly reduced abundance and diversity of meiofauna in subtidal areas close to the source of the spill (Montagna et al., 2013).
In contrast to information on subtidal meiofaunal populations, data on community composition and variation within intertidal meiofauna communities within the GOM are lacking. To the best of our knowledge, there have been no studies that examine variation and yearly composition within GOM intertidal meiofauna populations. In the wake of the DWHOS, Bik et al. (2012a) explored intertidal meiofauna communities along Dauphin Island and in Mobile Bay, Alabama, using high-throughput sequencing approaches. They showed a dramatic community shift from predominantly metazoan composition prior to the DWHOS to communities dominated by fungal taxa after the spill (Bik et al., 2012a). What was not known in this study was the degree of variability or changes in community composition over time (e.g., months, seasons, years).
Herein, we explore intertidal meiofaunal community variation along sites along Dauphin Island and in Mobile Bay, Alabama, over the course of a year. To facilitate comparison, we chose the same localities as Bik et al. (2012a). This study is the first within the GOM region to examine community composition of intertidal meiofauna populations over the course of a year and the first to incorporate Illumina platform technologies to assess meiofaunal diversity. Goals of this work included (1) assessing spatial and temporal variation within intertidal meiofaunal communities, (2) determining whether the large fungal composition seen in Bik et al. (2012a) was persistent, and (3) ascertaining whether the DWHOS has a lingering impact on intertidal communities. This study provides a fundamental baseline to examine similar intertidal meiofauna community impacts of future natural and anthropogenic disturbances in the GOM region.
Materials and Methods
Between July 2011 and July 2012, bi-monthly sediment samples were taken from the same five intertidal locations (Fig. 1 and Appendix Table Al) along Dauphin Island and in Mobile Bay, Alabama, as Bik et al. (2012a). Eight sediment cores 10 cm in depth were collected from each location at each sampling date using a PVC pipe with an interior diameter of 4 cm. All samples were covered by water at the time of collection. In the field, each core was divided into 0-3-cm and 3-10-cm sections. Six cores were immediately frozen on dry ice for genetic analysis and two cores were preserved in DMSO EDTA Salt Solution (DESS) (Yoder et al., 2006) for a nematode morphological study to be published later. Upon return to Auburn University, frozen samples were stored at -80[degrees]C until processed.
Meiofauna isolation and DNA extraction
For each sampling locality and time, three cores were haphazardly chosen for processing, and the two depth sections were processed separately. Depth fractions (either 0-3 cm or 3-10 cm) from the three cores were weighed and combined into one sample in an attempt to reduce the sample-specific heterogeneity typically found in meiofaunal communities. Meiofauna were isolated from sediment by decanting. Sediment was agitated by vigorously inverting the sample 10 times in a flask with a 32-ppt Instant Ocean solution (Blacksburg, VA) up to a total of 950 ml, the larger particles were allowed to settle (30-s maximum), and then the aqueous layer was carefully poured over a 45-[micro]m sieve. This decantation protocol was repeated 10 times per sample, except for Bayfront Park's (BP) 3-10-cm fractions, which were decanted only five times due to a large amount of small silt particles that hindered filtration. Material retained on the sieve was transferred to a 50-ml conical tube and stored at -80[degrees]C until DNA was extracted. Most of the material obtained with this protocol falls in the defined size range of meiofauna. However, some slightly larger organisms (<2 mm) might have been retained as well, and unpublished data (Brannock and Halanych) show that decanting preferentially retains metazoans compared to some other eukaryotic taxa. Herein, we use the term "meiofauna" to reflect the entire biotic sample retained on the sieve.
Total genomic DNA was extracted using the PowerSoil DNA Isolation Kit (MoBio Laboratories, Cat #12888). The bead solution was transferred from each PowerBead tube to a sterile 1.5-ml tube. Decanted material was mixed, and 1 ml was placed into the PowerBead tube and centrifuged at 15,700 X g for 1 min; most of the supernatant was discarded. The bead solution was returned to the PowerBead tubes, and the manufacturer's protocol was followed except the 10-min vortexing step was replaced with a 2-min Mini Beadbeater (BioSpec) step and DNA was eluted in a final volume of 50 [micro]l. DNA integrity was checked by gel electrophoresis, and DNA was stored at -20[degrees]C until sent for amplification and sequencing.
Sequencing and bioinformatics
Genomic DNA was sent to the Genomic Services Laboratory at Hudson Alpha Institute for Biotechnology (Huntsville, Alabama) for amplicon sequencing using the eukaryotic-specific hypervariable V9 region of the 18S small subunit ribosomal RNA (SSU rRNA) gene utilizing primers 1380F and 151 OR (Amaral-Zettler et al., 2009). Samples were amplified using 20 ng per sample per technical replicate for the input for PCR, except for where DNA concentration was extremely low. In such cases, the DNA volume sent for amplification was divided evenly between two PCR replicates. Maximum DNA volume per PCR was 18 [micro]l due to reaction limits. Each sample was amplified in duplicate. These technical replicates remained separated through operational taxonomic unit (OTU) clustering in the analysis pipeline. Paired-end 100-bp reads were obtained by sequencing dual barcoded amplicons on an Illumina HiSeq 2500 platform or paired-end 150-bp reads on an Illumina MiSeq. All samples had at least one technical replicate run on the HiSeq 2500, while only a subset of samples had one replicate run on the MiSeq. All samples from March 2012, May 2012, and July 2012 as well as Cadillac Ave 3-10-cm sample from July 2011 and BelleAir Blvd 3-10-cm sample from September 2011 had one technical replicate run on the MiSeq and one technical replicate run on the HiSeq 2500. All other samples had both technical replicates run on the HiSeq 2500. Previous studies have showed consistent results in amplicon analysis between the HiSeq and MiSeq platforms (Caporaso et al., 2012). Reads were demultiplexed by the Genomic Services Laboratory. Raw high-throughput sequencing reads were deposited in the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA) database (Accession Number SRP042047).
Raw forward and reverse reads were overlapped using PandaSeq ver. 2.5 (Masella et al, 2012). During the overlapping process, any reads not containing the 1380F and 1510R primers were discarded, and forward and reverse primers were trimmed after overlapping. In addition, sequences that contained uncalled bases (N) or were longer than 200 bp were also discarded. Remaining sequences were quality filtered utilizing the FASTX-Toolkit ver. 13.2 fastq_quality_trimmer script (Hannon Lab, 2009) with a quality score cutoff of 30. Chimeras were identified and filtered in QIIME ver. 1.8 (Caporaso et al., 2010b) using USEARCH ver. 6.1 (Edgar, 2010), and the remaining sequences were combined into one fasta file with properly formatted labels using the QIIME script add_qiime_labels.py. Sequences less than 75 bp were discarded.
Sequences were clustered into OTUs at 97% similarity using QIIME's pick_open_reference_otus.py workflow with the UCLUST (Edgar, 2010) clustering method disallowing singletons. Alternative clustering methods (e.g., USEARCH and UPARSE--Edgar, 2010 and 2013, respectively) were tested on another smaller dataset (PMB and KMH, unpubl. data) and did not differ greatly in the resulting taxonomic assignments; therefore, since the current dataset is too large to utilize the 32-bit ver. of USEARCH, clustering was conducted using UCLUST. The SILVA 111 (Quast et al., 2013) 99% clustered database was utilized as the reference set. The most abundant sequence was chosen as the reference sequence for the resulting OTU cluster. Taxonomy was assigned to representative sequences by using MegaBLAST (sequence identity [greater than or equal to] 90%, e-value 0.000001) to retrieve the top-scoring hit against the SILVA 111 99% clustered database. Alternative taxonomic assignment methods (e.g., RDP Classifier and UCLUST--Wang et al., 2007, and Edgar, 2010, respectively) did not result in OTU designations more specific than Eukaryota or Metazoa and were therefore not useful for these analyses. OTUs that did not result in a hit were classified as a "no hit." We restricted the sequence identity match to [greater than or equal to] 90% to be consistent with Bik et al. (2012a) methods as well as to be more conservative in our taxonomic assignments. Representative OTUs were aligned against the same reference database (SILVA 111, 99% clustered) using PYNAST (Caporaso et al., 2010a) (75-bp minimum length, 70% minimum ID). A minimum evolution tree was constructed with the resulting alignment using FASTTREE (Price et al., 2009). The resulting tree was utilized for alpha and beta-diversity measures (see the following section). Sequences classified as Bacteria and Archaea as well as those that failed to align with PYNAST were filtered from the resulting OTU table. Hereafter, the OTU table that excluded those sequences that failed to align with PYNAST and were classified as Bacteria or Archaea will be referred to as the final OTU table.
To determine whether technical replicates resulted in similar OTU composition, a hierarchical clustering analysis on the proportional OTU table for the 140 samples (70 samples with two PCR replicates each) was performed using the BiodiversityR package (Kindt and Coe, 2005) in R ver. 3.0.0 (R Development Core Team, 2008) utilizing an average linkage agglomerative clustering algorithm. The distance matrix was calculated using the Bray-Curtis dissimilarity metric. In the resulting dendrogram, samples with more similar OTU composition are clustered together and samples with more distinct OTU compositions cluster separately. As the dendrogram showed that technical replicates were more similar in OTU composition to each other than to any other sample (Appendix Fig. Al), sequences for technical replicates were combined and treated as one sample for the remainder of the analyses. Likewise, differences in replicate runs on the two platforms (Illumina HiSeq 2500 and MiSeq) were not observed. Hereafter, the data that include the combined replicates will be referred to as the complete dataset.
To allow direct comparison of diversity measures between samples, we normalized the number of sequences compared from each sample using a rarefaction analysis that subsampled the complete data set at 10 levels, from 10 to 300,010 sequences, in steps of 30,000 sequences. Ten replicates for each subsampling level were conducted, and average alpha-diversity measurements for each level were plotted against depth to determine a sufficient subsampling level. We subsampled the OTU table to 69,000 sequences in subsequent alpha and beta-diversity analyses. In order to be conservative in our community measurements, this number was chosen because it was close to the number of sequences (69,224) in our smallest sample that contained greater than 50,000 reads. For each subsample, three separate alpha-diversity measurements were calculated: Chaol, phylogenetic distance (PD), and Shannon diversity. Both Chaol and PD examine species richness, while Shannon diversity examines species evenness. Beta diversity was measured by calculating the weighted Unifrac distance (Lozupone and Knight, 2005), Bray-Curtis dissimilarity, and the binary Jaccard dissimilarity metrics in QIIME. The weighted Unifrac metric is a quantitative metric based on phylogenetic relationship of the OTUs that looks at the evenness of samples and takes into account the relative abundance of the OTUs. The Bray-Curtis metric measures species turnover on the basis of the raw counts of the OTUs, and the Jaccard metric is based on presence-absence of an OTU. Neither the Bray-Curtis nor the Jaccard metrics take into account species relatedness and therefore are not based on any phylogenetic information. Moreover, use of a presence-absence based approach allows assessment of OTU data irrespective of abundance or counts. In essence, the Jaccard metric allows one to more fully understand whether differences between samples represent changes in community composition by negating variation in abundances between OTUs due to differences in copy number of ribosomal genes. Jackknifed (n = 100) beta-diversity matrices were utilized in both principle coordinate analysis (PCoA) and unweighted pair group method with arithmetic mean (UPGMA) hierarchical clustering analysis.
Community composition was examined by summarizing taxonomic assignments with the summarize_taxa_through_ plots.py script in QIIME. This was performed on both the complete dataset and the subsampled OTU table normalized to 69,000 sequences. This study focused on the metazoan fraction of the eukaryotes, and for purposes of reporting, any metazoan phylum that had an average abundance less than 0.5% throughout the whole OTU table was classified as "Other Metazoa." The reasoning behind creating the "Other Metazoa" category was to reduce the number of very low abundance taxa that would otherwise render the bar graph unreadable. In addition, OTUs classified as fungi or stramenopiles, as well the "no hit" category, were reported separately when discussing the proportion of taxa. All other OTUs were grouped into an "Other Eukaryote" category that included unicellular eukaryotes. Phyla comprising the "Other Metazoa" and "Other Eukaryote" categories are listed in Appendix Table A2.
In all, 70 environmental samples were obtained representing five localities, seven time points, and two depth fractions. Five samples had extremely low initial DNA concentrations and therefore had less than 20 ng of input into the amplicon PCR. These included Shellfish Lab (SL) 0-3 cm from September 2011 and from November 2011, Ryan Court (RC) 0-3 cm from November 2011, BelleAir Blvd (BB) 0-3 cm from November 2011 and 3-10 cm from November 2011. All of these samples except BB 3-10 cm from November 2011 had greater than 100,000 sequences once the technical replicates were combined. BB 3-10 cm from November 2011 had only 1757 sequences after combining technical replicates.
The 70 environmental samples, each with two technical replicates, yielded a total of 24,381,423 demultiplexed Illumina reads in each paired-end direction. After preliminary overlapping (2,404,230 reads removed, <10% of original reads), quality (86 sequences removed, <0.0004% of overlapped sequences), chimera (5442 sequences removed, <0.03% of quality checked sequences), and length filtering (770 sequences removed, <0.004% of sequences after chimera checking), 90% (21,970,895 sequences) of original reads went into the OTU clustering process. Filtered sequences from the MiSeq ranged from 75 to 200 bp with both an average and median length of 138 bp. Sequences from the HiSeq 2500 ranged from 75 to 159 bp with an average of 138.5 bp and a median of 138 bp in length. After filtering, the final OTU table consisted of 98,765 OTUs across all samples. Removal of OTUs that failed to align with PYNAST (1953 OTUs, 1.9% of OTUs, containing 277,816 sequences, 1.3% of sequences) or that were classified as Bacteria or Archaea (656 OTUs, 0.6% of PYNAST filtered OTU table, containing 3164 sequences, 0.015% sequences after failed PYNAST alignment filtered) was limited. Of the 98,765 OTUs, only 7% (6931 OTUs containing 74,860 sequences, <0.4% of sequences in the final OTU table) were classified as "no hit." The "no hit" OTUs were blasted (blastn, e-value 0.000001) against the NCBI Genbank database, resulting in a majority (5971 OTUs, 86.1% of "no hit" OTUs and 6% of all OTUs in the final table) of reference sequences having no found match. The remaining 960 OTUs matched an rRNA entry in the database, with a majority (644 OTUs, 9.3% of "no hit" OTUs and <0.7% of all OTUs in the final table) matching an entry less than the initial 90% sequence identity cutoff. The remaining 316 "no hit" OTUs (4.6% of "no hit" OTUs and <0.4% of all OTUs in the final table) did match an NCBI entry with [greater than or equal to] 90% sequence identity. Discrepancies could be caused by the version of the SILVA database used (current QIIME compatible version is 111) or the fact that the SILVA database is curated and some NCBI entries were discarded during the curation process. The 70 samples ranged from 1,757 to 814,604 sequences, with a mean of 301,799 sequences.
Rarefaction analyses suggest that species richness (Choal and PD) were not saturated at a sampling depth of 300,010 sequences, although evenness (Shannon diversity) leveled off at greater than 30,000 sequences per sample (Appendix Fig. A2). Four samples (BP 3-10 cm on July 2011, November 2011, and July 2012 as well as BB 3-10 cm on November 2011) were below the subsampling threshold of 69,000 sequences and therefore were removed from further analyses.
Taxonomic assignments did not differ greatly between the OTUs based on the complete and subsampled datasets. These two datasets differed only by, at most, 0.3% for a given taxonomic grouping. Therefore, presentation of community taxonomic composition is based on the complete dataset. However, since the number of sequences in each OTU could vary depending on the sequence depth, the subsampled OTU table was used for all alpha- and beta-diversity measurements.
Communities examined were predominantly composed of metazoan taxa, and fungal OTUs were nonexistent or present in very low abundance (Fig. 2). A majority (94%) of communities had 72%-99% metazoan composition, while 6% (4 samples) contained 40%-63% metazoans. Annelids, arthropods, nematodes, and platyhelminths were the most abundant taxa throughout the samples (Fig. 2). Sample locations, however, differed in dominant taxa present. Locations were heterogeneous in amount of community variation observed over the year sampled. For example, RC and BB varied little in community composition in comparison to SL and BP (Fig. 2). Depth fractions (0-3 cm and 3-10 cm) differed in composition, and this difference varied temporally and spatially (Fig. 2). When comparing community composition between July 2011 and July 2012, some of the locations showed differences in taxa composition between years (Fig. 2), suggesting potential for yearly variation within the composition of these communities as well; however, additional data are required to test the yearly variation in these communities.
Alpha- and beta-diversity
For the three alpha-diversity measures, sample locations showed limited variation in alpha-diversity over the year sampled, except for SL (Appendix Fig. A3). Thus, even though community composition was variable in reference to the taxa present at some locations, the overall number of species present stayed relatively consistent.
Jackknifed UPGMA hierarchical clustering analysis (Fig. 3) and PCoA analysis (Fig. 4) revealed that samples showed separation by collection site. This clustering by site was apparent in all three beta-diversity matrices (Fig. 4), and Analysis of Similarity (ANOSIM) showed significant groupings of the beta-diversity of communities in reference to collection site (Table 1). However, the weighted unifrac showed a lower R-statistic in comparison to Jaccard or Bray-Curtis metrics (Table 1). When samples were grouped by location type, bay side (CA, BP, and BB) or Gulf side (RC and SL), there was a significant difference in the diversity of communities, and the finding was consistent for all three beta-diversity measures (Table 1). Samples did not appear to cluster by depth fraction or season (Figs. 3 and 4). ANOSIM results indicated that sample variation in relation to depth or season was significantly different from random (Table 1). Results for the Jaccard beta-diversity for season showed a significant P value (P = 0.011); however, the R-statistic was close to zero (R = 0.0712), indicating the group was not quite random. In addition, similar to the proportion of taxa breakdown, the two July sampling time points for a given location did not cluster together (Figs. 3 and 4), indicating potential yearly variation in community composition.
Illumina amplicon sequencing of the eukaryotic-specific hypervariable V9 SSU rRNA gene region proved useful in examining the community composition in intertidal meiofauna locations from Mobile Bay and Dauphin Island, Alabama. The high-throughput sequencing approach incorporated in this study showed variability in community composition over both spatial and temporal scales. Even though the number of species did not vary greatly throughout the data (alpha-diversity, Appendix Fig. A3), there was variability in which taxa were present. For example, during the September 2011 sampling of SL we noted an abundance of bryozoan biomass floating in the surf as well as washed up on the shore (Brannock, pers. obs.). In the community composition breakdown for SL 0-3 cm (Fig. 2) there was a bryozoan signature not found in other samples.
Overall results showed that annelids, arthropods, nematodes, and platyhelminthes were the most abundant taxa throughout the dataset; however, community composition showed spatial and temporal variation throughout the year (Fig. 2). RC and BB were heavily dominated by nematodes year round as well as in both depth fractions, and variation in community composition between sampling dates was low in comparison to SL (Fig. 2). Community composition was highly variable throughout the year at SL, especially for the 0-3-cm fraction, where from September 2011 through May 2012 community composition fluctuated greatly between several metazoan taxa (Fig. 2). On the other hand, the SL 3-10-cm fraction gradually progressed from a predominately annelid composition to a community composition of arthropod, nematode, and playhelminthes taxa during the same time frame. The CA 0-3-cm fraction appeared to have a fairly consistent nematode taxa proportion throughout the year, but varied in the amount of annelid, arthropod, brachiopod, and echinoderm taxa present. The 3-10-cm fraction showed variation in the proportions of major taxa present (annelids, arthropods, brachiopods, and echinoderms) and had a high nematode proportion except for July 2011, which was dominated by annelids. The BP 3-10-cm fraction consisted mainly of nematode and stramenopile taxa up to January 2012, at which time it drastically shifted to a predominately annelid community. After January 2012, annelid taxa decreased in this depth fraction, and the community returned to mainly nematode and stramenopile composition by July 2012. That same fluxuation in annelid taxa in January 2012 was present in the 0-3-cm fraction as well. Previous studies have shown meiofauna vary seasonally (Harris, 1972; Herman and Heip, 1983; Coull, 1985) and yearly (Harris, 1972; Herman and Heip, 1983; Coull, 1985; Pequegnat et al, 1990) in both species densities (Harris, 1972; Herman and Heip, 1983; Coull, 1985; Pequegnat et al., 1990) and vertical distribution (Harris, 1972); however, the degree of variation observed was species-dependent (Harris, 1972; Herman and Heip, 1983; Coull, 1985).
Community composition varied by location, as shown by samples clustered mainly by collection site (Figs. 3 and 4, Table 1). Likewise, samples collected from bay and gulf locations clustered separately from each other (Fig. 3 and Table 1). This suggests that location rather than season or depth potentially plays a larger role in community composition. These groupings were more apparent in Jaccard analyses, but were also seen in both the Bray-Curtis and weighted Unifrac analyses (Fig. 4 and Table 1), suggesting that locations may have different operational taxonomic units (OTUs) present in their communities; however, when comparing raw and relative abundances of taxa, the differences become less apparent. A corresponding morphological nematode study shows similar clustering of samples by location rather than season (Sharma, unpubl. data). Sediment composition, location characteristics, environmental variables, or surrounding vegetation could all be contributing to the differences in community composition observed. Differences in sediment size, type, and structure have been noted to influence meiofaunal community composition (Coull, 1988; Giere, 2009). Sediment analysis at these five locations using Bik et al.'s (2012a) samples as well as a March 2011 (H. Bik, University of Birmingham, UK; unpubl. data) time point shows differences in grain size, sorting ability, and amount of organic carbon (Williams, 2013). Locations and sample dates both show these differences (Williams, 2013). Likewise, Coull (1985) noted meiofauna in a muddy location had higher abundance variation in comparison to a more sandy location. More recently, Baguley et al. (2006) reported an increase in meiofauna abundance with increase in silt concentration. RC and SL are both more sandy coastal locations, while CA, BP, and BB are more protected bay locations that have higher silt content than RC and SL (Brannock, pers. obs.). Sediment composition and properties should be incorporated into future studies to examine the influence on community composition further. In addition, environmental factors such as temperature and salinity have also been reported as playing a role in the composition of meiofauna communities (Coull, 1985, 1999; Moens et al., 2013). Coull (1985) reported that salinity was correlated with meiofaunal abundance in the muddy location, while temperature was correlated with the sandy locations. During the current study, RC consistently had the highest salinity, whereas both BP and BB had relatively lower salinities in comparison to the other locations (Appendix Table Al). Temperature, on the other hand, showed less variation between the sample locations at a given time point (Appendix Table Al). Although environmental measurements were obtained only during sample collection, our limited number of time points for specific sampling localities precluded the ability to meaningfully correlate salinity and temperature data to variations observed in OTU data. To have a better understanding of the influences of these variables on the community composition, more regular measurements need to be obtained to grasp the daily and monthly fluctuations.
For Gulf of Mexico (GOM) meiofaunal localities sampled here, by July 2011 communities were similar to pre-spill communities (fig. 1 in Bik et al., 2012a) in that metazoan taxa dominated the composition (Fig. 2). In contrast, samples in September 2010 were dominated by fungal species (Bik et al., 2012a). By July 2011, fungal OTUs were not detectable or were present in extremely low abundances (Fig. 2), suggesting that the drastic shift in community composition that Bik et al. (2012a) reported as a potential impact of the Deepwater Horizon oil spill (DWHOS) was an acute phenomenon. Even though fungal impact on these sites might have been short term, communities could be experiencing prolonged impacts from the DWHOS disturbance. When combining the two depth fractions (0-3 cm and 3-10 cm) in this study for the May 2012 time point (Fig. 5) and comparing the taxa composition to those reported in Bik et al. (2012a) for May 2010, community compositions differ between the two years. One of the most notable differences is that the RC May 2012 community composition of predominantly nematode and arthropod taxa (Fig. 5) resembled the post-spill community in Bik et al. (2012a, Fig. 1) more closely than it resembled the May 2010 composition of predominately annelid taxa.
There were methodological differences between Bik et al. (2012a) and the current study that could also potentially explain the observed variation in community compositions between the two studies. First, Bik et al. (2012a) used a smaller amount (200 jul) of decanted material in their DNA extractions and used a different extraction kit. We attempted to use 200 p\ of decanted material in our extractions at first with little success in DNA yield. Therefore, we increased the amount of initial extraction material. However, when using a larger amount of starting material one would not expect to find taxa absent among the fungi. In addition, the 18S gene region targeted differed between the two studies. Bik et al. (2012a) used primer pairs that targeted a larger fragment (--450 bp) toward the beginning of the 18S gene, while we targeted a smaller fragment (87-187 bp, Amaral-Zettler et al., 2009) toward the end of the 18S gene. Primers used in the current study were chosen because at the time of the experiment the Illumina platform was unable to sequence through the full fragment used in Bik et al. (2012a).
Determining whether the Deepwater Horizon oil spill (DWHOS) may have a long-term impact on meiofaunal communities is difficult because knowledge of seasonality and variability in these communities is still limited. The current study did illustrate seasonal and yearly variation in community composition at the sample locations. However, whether this variation is normal or indicative of a community re-stabilizing after the DWHOS disturbance is unknown. Continued sampling over a longer timescale as well as the inclusion of environmental parameters is necessary to comprehend the drivers in variation observed in these communities.
Table A1 The Mobile Bay and Dauphin Island, Alabama, sample locations with the corresponding abbreviations and GPS coordinates as well as salinity and water temperature measurements taken at the time of sampling (image of the sample locations can be found in Fig. 1) Ryan Cadillac Shellfish Court Ave Lab Abbreviation RC CA SL Latitude 30.2500N 30.2532N 30.2465N Longitude 88.1462W 88.1363W 88.0784W Salinity (ppt) Jul 2011 29.0 20.0 22.0 Sep 2011 25.0 20.0 20.0 Nov 2011 30.0 28.0 28.0 Jan 2012 24.0 19.0 10.0 Mar 2012 32.0 24.5 27.0 May 2012 27.0 17.0 23.0 Jul 2012 29.0 24.0 23.0 Water Jul 2011 N/A N/A N/A Temperature Sep 2011 30.8 30.8 28.3 ([degrees]C) Nov 2011 N/A N/A N/A Jan 2012 18.6 20.4 19.8 Mar 2012 22.2 23.2 22.7 May 2012 29.3 27.4 29.2 Jul 2012 31.8 32.5 33.1 Bayfront BelleAir Park Blvd Abbreviation BP BB Latitude 30.3538N 30.5080N Longitude 88.1179W 88.1016W Salinity (ppt) Jul 2011 17.0 10.0 Sep 2011 10.0 7.0 Nov 2011 20.0 15.0 Jan 2012 8.0 5.0 Mar 2012 5.0 4.0 May 2012 14.0 11.0 Jul 2012 18.0 11.0 Water Jul 2011 N/A N/A Temperature Sep 2011 29.4 26.0 ([degrees]C) Nov 2011 N/A N/A Jan 2012 22.6 21.4 Mar 2012 23.4 23.9 May 2012 29.5 28.2 Jul 2012 35.0 33.7 N/A, not available. Table A2 List of phyla comprising either the "Other Metazoa" or the "Other Eukaryotes" groupings in Figs. 2 and 5 Other Metazoa Phyla Other Eukaryote Phyla Chaetognatha Ciliophora Ctenophora Other Alvelata Entoprocta Amoebozoa Hemichordata Cryptophyta Loricifera Loukozoa Gnathifera Glaucophyta Onochophora Choanozoa Porifera Picozoa Rotifera Foraminifera Xenoturbeiidea Rhodophyta Acanthocephala Myzozoa Chordata Apusozoa Cycliophora Heliozoa Gnathostomulida Percolozoa Kinohycha Euglenozoa Mesozoa Haptophyta Nematomorpha Metamunada Placozoa Cercozoa Priapulida Radiozoa Sipuncula Uncultured unicellular Uncultured Metazoa Eukaryote
We acknowledge Dr. Holly M. Bik, Stephanie K. Hoffman, Dr. Molli M. Newman, and Dr. Scott R. Santos for all their support and assistance throughout this project. We also thank David R. Branson, Matthew P. Galaska, Dr. Alexis M. Janosik, Kiley W. Seitz, and Amanda O. Shaver for their assistance with sample collection. This research was funded by Year 1 Gulf of Mexico Research Initiative (GoMRI) awarded by Marine Environmental Science Consortium (MESC). This is Molette Biology Laboratory contribution 26 and Auburn University Marine Biology Program contribution 119.
Received 13 March 2014; accepted 17 June 2014.
Aller, R. C., and J. Y. Alter. 1992. Meiofauna and solute transport in marine muds. Limnol. Oceanogr. 37: 1018-1033.
Amaral-Zettler, L. A., E. A. McCliment, H. W. Ducklow, and S. M. Huse. 2009. A method for studying protistan diversity using massively parallel sequencing of V9 hypervariable regions of small-subunit ribosomal RNA genes. PLoS One 4: e6372.
Baguley, J. G., P. A. Montagna, L. J. Hyde, R. D. Kalke, and G. T. Rowe. 2006. Metazoan meiofauna abundance in relation to environmental variables in the Northern Gulf of Mexico deep sea. Deep Sea Res. Part I 53: 1344-1362.
Baguley, J. G., P. A. Montagna, L. J. Hyde, and G. T. Rowe. 2008. Metazoan meiofauna biomass, grazing, and weight-dependent respiration in the Northern Gulf of Mexico deep sea. Deep Sea Res. Part II 55: 2607-2616.
Bik, H. M., K. M. Halanych, J. Sharma, and W. K. Thomas. 2012a. Dramatic shifts in benthic microbial eukaryote communities following the Deepwater Horizon oil spill. PloS One 7: e38550.
Bik, H. M., D. L. Porazinska, S. Creer, J. G. Caporaso, R. Knight, and W. K. Thomas. 2012b. Sequencing our way towards understanding global eukaryotic biodiversity. Trends Ecol. Evol. 27: 233-243.
Bik, H. M., W. Sung, P. De Ley, J. G. Baldwin, J. Sharma, A. Rocha-Olivares, and W. K. Thomas. 2012c. Metagenetic community analysis of microbial eukaryotes illuminates biogeographic patterns in deep-sea and shallow water sediments. Mol. Ecol. 21: 1048-1059.
Caporaso, J. G., K. Bittinger, F. D. Bushman, T. Z. DeSantis, G. L. Andersen, and R. Knight. 2010a. PyNAST: A flexible tool for aligning sequences to a template alignment. Bioinformatics 26: 266267.
Caporaso, J. G., J. Kuczynski, J. Stombaugh, K. Bittinger, F. D. Bushman, E. K. Costello, N. Fierer, A. G. Pena, J. K. Goodrich, J. I. Gordon et al. 2010b. QIIME allows analysis of high-throughput community sequencing data. Nat. Methods 7: 335-336.
Caporaso, J. G., C. L. Lauber, W. A. Walters, D. Berg-Lyons, J. Huntley, N. Fierer, S. M. Owens, J. Betley, L. Fraser, M. Bauer et al. 2012. Ultra-high-throughput microbial community analysis on the Illumina HiSeq and MiSeq platforms. ISME J. 6: 1621-1624.
Coull, B. C. 1985. Long-term variability of estuarine meiobenthos: an 11 year study. Mar. Ecol. Prog. Ser. 24: 205-218.
Coull, B. C. 1988. Ecology of the marine meiofauna. Pp. 18-38 in Introduction to the Study of Meiofauna, R. P. Higgins and H. Thiel, eds. Smithsonian Insitution Press, Washington, DC.
Coull, B. C. 1999. Role of meiofauna in estuarine soft-bottom habitats. Aust. J. Ecol. 24: 327-343.
Creer, S., V. G. Fonseca, D. L. Porazinska, R. M. Giblin-Davis, W. Sung, D. M. Power, M. Packer, G. R. Carvalho, M. L. Blaxter, P. J. D. Lambshead, and W. K. Thomas. 2010. Ultrasequencing of the meiofaunal biosphere: practice, pitfalls and promises. Mol. Ecol. 19: 4-20.
Edgar, R. C. 2010. Search and clustering orders of magnitude faster than BLAST. Bioinformatics 26: 2460-2461.
Edgar, R. C. 2013. UPARSE: Highly accurate OTU sequences from microbial amplicon reads. Nat. Methods 10: 996-998.
Escobar, E., M. Lopez, L. A. Soto, and M. Signoret. 1997. Density and biomass of the meiofauna of the upper continental slope in two regions of the Gulf of Mexico. Cienc. Mar. 23: 463-489.
Escobar-Briones, E. G., and L. A. Soto. 1997. Continental shelf benthic biomass in the Western Gulf of Mexico. Cont. Shelf Res. 17: 585-604.
Fonseca, V. G., G. R. Carvalho, W. Sung, H. F. Johnson, D. M. Power, S. P. Neill, M. Packer, M. L. Blaxter, P. J. D. Lambshead, W. K. Thomas, and S. Creer. 2010. Second-generation environmental sequencing unmasks marine metazoan biodiversity. Nat. Commun. 1: 98.
Giere, O. 2009. Meiobenthology: The Microscopic Motile Fauna of Aquatic Sediments. 2nd ed. Springer-Veriag, Berlin
Hannon Lab. 2009. FAST-X Toolkit: FASTQ/A short-reads pre-processing tools. [Online]. Cold Spring Harbor Laboratory, Cold Spring Harbor, NY. Available http://hannonlab.cshl.edu/fastx_toolkit/ [2014 September 17].
Harris, R. P. 1972. Seasonal changes in the meiofauna population of an intertidal sand beach. J. Mar. Biol. Assoc. UK 52: 389-403.
Herman, P. M. J., and C. Heip. 1983. Long-term dynamics of meiobenthic populations. Proceedings 17th European Marine Biology Symposium, Brest, France, 27 September-1 October 1982. Oceanol. Acta 1983: 109-112.
Kindt, R., and R. Coe. 2005. Tree Diversity Analysis: A Manual and Software for Common Statistical Methods for Ecological and Biodiversity Studies. World Agroforestry Centre (ICRAF), Nairobi.
Landers, S. C., F. A. Romano III, P. M. Stewart, and S. Rampoop. 2012. A multi-year survey of meiofaunal abundance from the Northern Gulf of Mexico continental shelf and slope. Gulf Mex. Sci. 30: 20-29.
Lozupone, C., and R. Knight. 2005. UniFrac: a new phylogenetic method for comparing microbial communities. Appl. Environ. Microbiol. 71: 8228-8235.
Mare, M. F. 1942. A study of a marine benthic community with special reference to the micro-organisms. J. Mar. Biol. Assoc. UK 25: 517-554.
Masella, A. P., A. K. Bartram, J. M. Truszkowski, I). G. Brown, and J. D. Neufeld. 2012. PANDAseq: Paired-end assembler for Illumina sequences. BMC Bioinformatics 13: 31.
Moens, T., U. Bracckman, S. Derycke, G. Fonseca, F. Gallucci, R. Gingold, K. Guilini, J. Ingels, D. 1. Ledue, J. Vanaverbeke et al. 2013. Ecology of free-living marine nematodes. Pp. 109-152 in Handbook of Zoology, Vol. 3, A. Schmidt-Rheas, ed. De Gruyter, Berlin.
Montagna, P. A., and D. J. Harper. 1996. Benthic infaunal long-term response to offshore production platforms in the Gulf of Mexico. Can. J. Fish. Aquat. Sci. 53: 2567-2588.
Montagna, P. A., and R. D. Kalke. 1992. The effect of freshwater inflow on meiofaunal and macrofaunal populations in the Guadalupe and Nueces Estuaries, Texas. Estuaries 15: 307-326.
Montagna, P. A., R. D. Kalke, and C. Ritter. 2002. Effect of restored freshwater inflow on macrofauna and meiofauna in upper Rincon Bayou, Texas, USA. Estuaries 25: 1436-1447.
Montagna, P. A., J. G. Baguley, C. Cooksey, I. Hartwell, L. J. Hyde, J. L. Hyland, R. D. Kalke, L. M. Kracker, M. Reuscher, and A. C. E. Rhodes. 2013. Deep-sea benthic footprint of the Deepwater Horizon blowout. PloS One 8: e70540.
Pequegnat, W. E., B. J. Gallaway, and L. H. Pequegnat. 1990. Aspects of the ecology of the deep-water fauna of the Gulf of Mexico. Am. Zool. 30: 45-64.
Price, M. N., P. S. Dehal, and A. P. Arkin. 2009. FastTree: computing large minimum evolution trees with profiles instead of a distance matrix. Mol. Biol. Evol. 26: 1641-1650.
Quast, C., E. Pruesse, P. Yilmaz, J. Gerken, T. Schweer, P. Yarza, J. Peplies, and F. O. Glockner. 2013. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41: D590-596.
R Development Core Team. 2008. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. ISBN 3-900051-07-0. [Online]. Available at http://www.Rproject.org.
Snelgrove, P. V. R. 1997. The importance of marine sediment biodiversity in ecosystem processes. Ambio 26: 578-583.
Wang, Q., G. M. Garrity, J. M. Tiedje, and J. R. Cole. 2007. Naive bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl. Environ. Microbiol. 73: 5261-5267.
Williams, S. C. 2013. Characterization of sediment around Dauphin Island, Alabama for a biological context of ecosystem recovery related to the BP Horizon oil spill. Masters thesis, The University of Texas at San Antonio.
Yoder, M., I. Tandingan De Ley, I. W. King, M. Mundo-Ocampo, J. Mann, M. Blaxter, L. Poiras, and P. De Ley. 2006. DESS: a versatile solution for preserving morphology and extractable DNA of nematodes. Nematology 8: 367-376.
PAMELA M. BRANNOCK (1) *, DAMIEN S. WAITS (1), JYOTSNA SHARMA (2), AND KENNETH M. HALANYCH (1)
(1) Department of Biological Sciences, Auburn University, Auburn, Alabama 36849; and (2) Department of Biology, University of Texas at San Antonio, Texas 78249
* To whom correspondence should be addressed. E-mail: pmb0010@ auburn.edu
Table 1 Analysis of Similarity (ANOSIM): resulting ANOSIM R-statistic * and P value for all three beta-diversity distance matrices (f) Bray-Curtis Jaccard R-statistic P R-statistic p Location 0.6422 0.001 0.7219 0.001 Type Depth 0.0133 0.204 0.0241 0.102 Season 0.0212 0.185 0.0712 0.011 Site 0.8129 0.001 0.8547 0.001 Weighted Unifrac R-statistic P Location 0.2777 0.001 Type Depth -0.0142 0.846 Season 0.0050 0.406 Site 0.4817 0.001 * An R-statistic of 0 represents complete random grouping, a value of +1 indicates samples are most similar within the same group, and a value of -1 indicates samples most similar outside of the group. ([dagger]) For each distance matrix, four separate comparisons were made to determine whether significant differences were seen between groups.
|Printer friendly Cite/link Email Feedback|
|Author:||Brannock, Pamela M.; Waits, Damien S.; Sharma, Jyotsna; Halanych, Kenneth M.|
|Publication:||The Biological Bulletin|
|Date:||Oct 1, 2014|
|Previous Article:||Detection and removal of PCR duplicates in population genomic ddRAD studies by addition of a Degenerate Base Region (DBR) in sequencing adapters.|
|Next Article:||Ciliate diversity, community structure, and novel taxa in lakes of the McMurdo Dry Valleys, Antarctica.|