Meiofaunal diversity at a seamount in the pacific ocean: A comprehensive study using environmental DNA and RNA

Deep–Sea Research I xxx (xxxx) xxx

Contents lists available at ScienceDirect

Deep-Sea Research Part I journal homepage: http://www.elsevier.com/locate/dsri

Meiofaunal diversity at a seamount in the pacific ocean: A comprehensive study using environmental DNA and RNA Tomo Kitahashi a, *, Sachie Sugime a, Kentaro Inomata a, b, Miyuki Nishijima a, b, d, Shogo Kato c, Hiroyuki Yamamoto a a

Japan Agency for Marine-Earth Science and Technology (JAMSTEC), 2-15 Natsushima-cho, Yokosuka, Kanagawa, 237-0061, Japan TechnoSuruga Laboratory Co., Ltd, 330 Nagasaki, Shimizu-ku, Shizuoka, Shizuoka, 424-0065, Japan Japan Oil Gas, And Metals National Corporation (JOGMEC), Toranomon Twin Building 2-10-1 Toranomon, Minato-ku, Tokyo, 105-0001, Japan d National Institute of Advanced Industrial Science and Technology (AIST), AIST Tsukuba Central 7, 1-1-1 Higashi, Tsukuba, Ibaraki, 305-8567, Japan b c

A R T I C L E I N F O

A B S T R A C T

Keywords: Seamount Meiofauna Metabarcoding analysis eDNA eRNA

In deep-sea environments, which have low densities of large benthic animals, meiofauna are a useful indicator of biodiversity. However, in low-latitudes and in low-productivity areas, meiofaunal density is low, but is higher than that of megafaunal or macrofaunal density. Therefore, it is difficult to collect a sufficient number of specimens for statistical analysis. In addition, because faunal classification has been based on conventional microscopic observations, the use of meiofauna to estimate biodiversity in deep-sea environments is timeconsuming. However, metabarcoding analyses focused on environmental DNA (eDNA) or RNA (eRNA) have recently been used to examine deep-sea eukaryotic diversity and communities. Here, we examined meiofaunal assemblages using microscopic, eDNA-, and eRNA-based methods at Xufu Guyot (JA06 Seamount), off the southeastern coast of Minami-Torishima Island in the North Pacific Ocean. Microscopic analysis failed to detect a significant difference in diversity or community structure between the seamount terrace and base. This was likely because of the low abundance of meiofauna, which was caused by the low surface productivity at the study area. However, eDNA/eRNA-based metabarcoding analyses revealed spatial variations in diversity and community structures within a single seamount. Therefore, metabarcoding analysis might be useful to elucidate meiofaunal assemblages in areas with low productivity and low faunal density.

1. Introduction Seamounts are oases of productivity and hotspots of biodiversity, making them important environments in deep-sea ecosystems (Rowden et al., 2010; George et al., 2018). Megafaunal and macrofaunal density and diversity are higher at seamount summits than at seamount slopes or comparable areas (Rowden et al., 2010; Sautya et al., 2011). There are 232 seamounts that have been surveyed for biological studies; however, meiofauna had only been surveyed at eight seamounts (George, 2013), although they are more abundant than macro and megafauna (Wei et al., 2010) and important components of deep-sea environments (e.g. Giere, 2009). There are two factors causing the limited knowledge on seamount meiofaunal assemblage. The first factor is difficulties in quantitative sampling and identification. Varieties of bottom structures (slope, even, and uneven) and bottom types (sediment, sand, and rock) make it difficult to quantitatively collect meiofaunal samples (George,

2013). The other one is research difficulties in identifying the complete species composition in a community, which are caused by the complexity of morphological characteristics. The solution is DNA signature for identifying taxon or phylogenetic position. Environmental DNA (eDNA) and environmental RNA (eRNA) can be used to analyze the biodiversity of the phylogenetic assemblages of prokaryotes and small eukaryotes such as meiofauna and foraminifera, although many novel DNA signature emerges without taxonomical nomenclature. The reverse-taxonomy approach using genetic analysis enables straightforward classification of individual genotypes instead of taxon names (Markmann and Tautz, 2005). Recent advances in sequencing technology have afforded metagenomic methods (also referred to as “metabarcoding”) based on high-throughput sequencing, which has been commonly used for evaluating deep-sea eukaryotic di­ versity and communities (Pawlowski et al., 2011; Bik et al., 2012; Guardiola et al., 2016).

* Corresponding author. E-mail address: [email protected] (T. Kitahashi). https://doi.org/10.1016/j.dsr.2020.103253 Received 1 August 2019; Received in revised form 14 February 2020; Accepted 17 February 2020 Available online 20 February 2020 0967-0637/© 2020 Elsevier Ltd. All rights reserved.

Please cite this article as: Tomo Kitahashi, Deep–Sea Research I, https://doi.org/10.1016/j.dsr.2020.103253

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

The emerging issue is the discrepancy between the results from eDNA and eRNA (Pochon et al., 2017). DNA extracted from environ­ mental samples contains detailed information on the structure and functions of in situ communities, including prokaryotes and small eu­ karyotes, although physiologically active, dormant, and dead cells are simultaneously present in such samples. DNA, which has been preserved for a long time in deep-sea sediments, can maintain its genetic infor­ mation (Dell’Anno and Danovaro, 2005; Corinaldesi et al., 2011). Therefore, eDNA approaches does not exclude data from the dead cells remained in the environment. In contrast, eRNA approaches targeting ribosomal RNA can distinguish between dead and live cells, as RNA decomposes quickly after cells die (Yamamoto et al., 1996; Takishita et al., 2010; Lejzerowicz et al., 2013). Thus, the common expectation that RNA is the most ephemeral molecule in environment, may change. Recent evidence suggest that RNA may be abundantly excreted by or­ ganisms and is sufficiently persistent in the environment to reconstruct community composition and gene expression (Cristescu, 2019). This assumption, however, remains largely untested in environmental studies, because comparative analysis with the conventional identifi­ cation method using microscopic characterization of meiofaunal morphology requires time-consuming examination. In this study, we examined meiofaunal community compositions using co-extracted eDNA/eRNA collected from seamount habitats, and compared the re­ sults with those obtained using the conventional method. The study was conducted at Xufu Guyot (JA06 Seamount) located in the high seas off Minami-Torishima Island in the Western Pacific Ocean (Fig. 1). We anticipated that the seamount formation from bottom to

summit would provide habitats for various meiofauna, and that it may be a suitable sampling site for this comparative study of meiofauna communities. In this study, we addressed the following two questions: 1) Are there any differences in meiofaunal assemblages in JA06 Seamount? 2) Is there any discrepancy between the results obtained using the conventional microscopic method and those obtained using eDNAand eRNA-based metabarcoding methods? 2. Materials and methods 2.1. Sampling and sample processing The R/V Hakurei underwent a cruise conducted by the Japan Oil, Gas and Metals National Corporation (JOGMEC) in May 2016. Sampling and observation were carried out at four stations on the base (B02, B02–2, B03, and B04) and two stations on the terrace (T02 and T02-2) of JA06 Seamount (Fig. 1b, Table 1). Sediment core samples were collected using a multiple corer (inner diameter: 7.4 cm) that could simultaneously collect up to eight sediment cores (Barnett et al., 1984). At B02, the collected sediments comprised diatomaceous ooze and overflew from the core tubes upon sample retrieval. Therefore, we collected additional sediment samples nearby this station (B02-2). Except for T02-2, one core was used for microscopic analysis, and the second core was used for metabarcoding analysis. At T02-2, three cores were used for microscopy, and one core was used for metabarcoding analysis. Samples for

Fig. 1. Locations of the study sampling stations. 2

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

with Ex Taq polymerase (TaKaRa Bio). The PCR conditions were as follows: initial denaturation for 1 min at 96 � C; 15 cycles for 30 s at 96 � C, 45 s at 65 � C, and 1 min at 72 � C; and final extension for 7 min at 72 � C. After purification with magnet beads, indexed PCR products were quantified using a Qubit dsDNA HS assay kit and a Qubit 2.0 Fluorom­ eter (Invitrogen). All purified index PCR products were diluted to the same concentration and pooled to create one metabarcoding library. Paired-end sequencing (600 bp) was subsequently performed using a MiSeq sequencer (Illumina) using a MiSeq Reagent Kit V3 (Illumina) according to the manufacturer’s protocol.

Table 1 Locations of the sampling stations. Station

Date

Latitude

Longitude

Depth (m)

Setting

B02 B02-2 B03 B04 T01 T02-2

28-May-2016 30-May-2016 27-May-2016 24-May-2016 25-May-2016 25-May-2016

19 23 19� 230 19� 230 19� 450 19� 320 19� 320

157 43 157� 430 158� 100 157� 550 157� 510 157� 000

4508 4417 4457 4516 1327 1312

Base Base Base Base Terrace Terrace

�

0

N N N N N N

�

0

E E E E E E

microscopic analysis were sliced horizontally into five 1-cm layers (i.e., 0–1, 1–2, 2–3, 3–4, and 4–5 cm), fixed, and then preserved individually in 5% buffered seawater formalin. Samples for metabarcoding analysis were collected from the 0–1 cm layer and preserved at 80 � C.

2.4. Sequencing analyses For raw forward and reverse reads, the primer region (i.e., SSU_F04/ SSU_R22mod) was trimmed using cutadapt (Martin, 2011) and the trimmed reads were overlapped with PEAR (Zhang et al., 2014). Sub­ sequent analyses were performed using QIIME version 1.9.1 (Caporaso et al., 2010b). Low-quality reads (QV < 20) were filtered using split_li­ braries_fastq.py, and chimeric reads were identified using identi­ fy_chimeric_seqs.py with the usearch61 method and removed using filter_fasta.py. Reads were clustered into operational taxonomic units (OTUs) at 97% similarity using pick_open_reference_otus.py against the reference database (SILVA 128, 97% clustered). In this step, the reads were clustered against the reference database, and the reads that did not match the database were subsequently clustered based on the similarity, and OTUs that contained less than 10 reads were removed using pick_­ open_reference_otus.py. We used SILVA 128 (97% clustered) as the reference database, which is the quality-controlled database of small subunit rRNA gene for bacteria, archaea, and eukaryotes (; Quast et al., 2013) and is frequently used in the metabarcoding studies of meiofauna (e.g. Bik et al., 2012; Sinniger et al., 2016). The representative reads of each OTU were aligned using align_seqs.py against the same reference database using PYNAST (Caporaso et al., 2010a), and the reads that could not be aligned were removed using filter_fasta.py. Next, these se­ quences were taxonomically assigned by assign_taxonomy.py according to the BLAST method (e-value: 1.0 6) against the data from SILVA 128 (97% clustered).

2.2. Microscopic analysis For microscopic analysis, Rose Bengal (final concentration: 0.05 g/L) was added to the preserved samples to stain the meiofaunal specimens. Next, the stained sediment samples were treated as previously described (Danovaro, 2010). The samples were passed through a 0.5-mm mesh sieve and retained on a 38-μm mesh sieve. The fraction remaining on the latter sieve was resuspended and centrifuged at 800 g three times in colloidal silica (Ludox HS-40, Sigma-Aldrich) for 10 min each time. The supernatants were transferred to flat-bottomed Petri dishes, and Rose Bengal-stained organisms were collected using an Irwin loop and iden­ tified up to high taxonomic taxa. 2.3. Metabarcoding analysis DNA was extracted from sediment samples by using a PowerMax Soil DNA Isolation kit (MO BIO Laboratories) according to the manufacture’s instruction with minor modifications. The V1 to V2 regions of the nu­ clear small subunit ribosomal RNA gene (18S rRNA gene) were ampli­ fied from the extracted DNA using the following primers: SSU_F04, 50 GCTTGTCTCAAAGATTAAGCC-30 ; SSU_R22mod, 50 -CCTGCTGCCT TCCTTRGA-30 (Sinniger et al., 2016). Both primers had an overhang adaptor sequence (forward primer: ACACTCTTTCCCTA­ CACGACGCTCTTCCGATCTGCTTGTCTCAAAGATTAAGCC; reverse primer: GTGACTGGAGTTCAGACGTGTGCTCTTCCGATCTCCTGCTGCC TTCCTTRGA) for index polymerase chain reaction (PCR) at each of their 50 end. PCR amplification was performed using 1 μL of genomic DNA template (1 ng/μL) in 20 μL reaction volume with Ex Taq polymerase (TaKaRa Bio). The PCR conditions were as follows: initial denaturation for 1 min at 96 � C; 25 cycles for 25 s at 96 � C, 45 s at 55 � C, and 1 min at 72 � C; and final extension for 7 min at 72 � C. RNA was extracted from sediment samples by using a PowerMax Soil RNA Isolation kit (MO BIO Laboratories). Extracted RNA was purified with RNA Clean & Concentrator-5 with DNase I Kit (ZYMO Research) to remove co-extracted DNA, and subsequently checked for any residual DNA by PCR and 2% agarose gel electrophoresis. Purified RNA was translated into complementary DNA (cDNA) by using reverse tran­ scriptase (ReverTra Ace, TOYOBO) and the reverse primer (SSU_R22­ mod) for 18S rRNA gene according to the manufacturer’s protocol. The synthesized cDNA was used as a template for amplification of the 18S rRNA gene with the same primers and procedure used for DNA samples. Three amplifications were performed in parallel for each station and for both DNA and RNA samples, and negative controls (i.e., ultrapure water only) were included in all amplification reactions. Triplicates of PCR products from DNA and RNA samples were visu­ alized in 2% agarose gel, and the amplified 500-bp fragment was puri­ fied using magnet beads (Agencourt AMPure XP Kit, Beckman Coulter). The purified PCR amplicons were attached to the sequencing adapters with sample index sequences of the MiSeq, and index PCR was per­ formed using 1 μL of PCR product (1 ng/μL) in 20 μL reaction volume

2.5. Diversity analyses For microscopic data, the following diversity indices were evaluated: number of taxa, Shannon diversity (H0 ), and Pielou’s evenness (J0 ) based on density data (ind./10 cm2). For metabarcoding data, the following diversity indices were evaluated: number of OTUs assigned to Metazoa, H0 , and J0 based on the present/absent data for each OTU. In addition, two diversity indices were evaluated for the metabarcoding data on the basis of the numbers of OTUs: 1) nonparametric richness estimator Chao 1 (Chao et al., 2005), which was used to estimate the asymptotic number of OTUs; and 2) the expected number of OTUs for a theoretical sample of 5000 reads (“Expected”), which was estimated using the rarefaction method (Hurlbert, 1971). 2.6. Statistical analysis Student’s t-test was used to determine the significance of differences in the meiofaunal density (microscopic data) and diversity indices (microscopic and metabarcoding data) between the base and terrace of the studied seamount. Multivariate analyses were performed to test the significance of differences in community structure for both microscopic and metabarcoding data. Principal co-ordinates (PCO) analysis was performed to visualize differences in community structure in meiofaunal analysis. Permutational multivariate analysis of variance (PERMA­ NOVA) was conducted to elucidate differences in community structure. Similarity profiles (SIMPROF) were determined to test the null hy­ pothesis that the community structures at each station did not differ and to identify statistically indistinguishable stations. 3

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

All statistical analyses were conducted with R studio (RStudio Team, 2016) using the vegan package (Oksanen et al., 2018).

Table 3 Meiofaunal diversity indices at the high taxonomic level in the 0–5 cm sediment layer at each station. Indices in the 0–1 cm sediment layer were presented in the parentheses.

3. Results 3.1. Microscopic analysis

B02 B02-2 B03 B04 T01 T02-2

Meiofaunal composition and diversity indices at each station are shown in Tables 2 and 3, respectively. Nematoda (including Desmo­ scolecid nematodes) was the most dominant taxon at all stations, and Copepoda was the second-most abundant taxon in all stations, except for B02 and B02-2 (Table 2). The meiofaunal density in terrace stations tended to be higher than that in base stations (Table 2), although the differences were not significant (t-test, p > 0.05). There were no sig­ nificant differences in diversity indices between base and terrace sta­ tions for both the 0–1 cm and 0–5 cm sediment layers (t-test, p > 0.05, Table 3). At all stations, except for B02 and B03, the meiofauna exhibited a typical vertical distribution pattern, i.e., meiofauna were aggregated at the surface layer of the sediment (Fig. 2). At B02, the maximum meio­ faunal density was observed at the subsurface layer (3–4 cm). Because the sediments samples at B02 overflowed from the core tubes upon sample retrieval, the vertical profile was disrupted at this station. At B03, copepods were dominant at the surface layer, with the maximum density at the subsurface layer. A PCO plot based on the assemblage composition at the high taxo­ nomic level is shown in Fig. 3. There was no apparent difference in assemblage composition between base and terrace stations, and there was no significant difference between settings for both the 0–1 cm and 0–5 cm sediment layers (PERMANOVA, p > 0.05).

Num of taxa

H0

4 (1) 7 (7) 3 (3) 6 (6) 7 (7) 8 (7)

1.0 1.3 0.7 1.1 1.1 0.6

J0 (0.0) (1.4) (1.0) (1.1) (1.2) (0.7)

0.71 ( ) 0.68 (0.71) 0.64 (0.86) 0.59 (0.63) 0.55 (0.63) 0.30 (0.35)

For both the DNA and RNA datasets, the numbers of reads and OTUs assigned to Cnidaria, Arthropoda, Nematoda, and Annelida were greater than those assigned to other phyla (Table 5). Among Cnidaria, Hydrozoa and Siphonophorae were dominant (Fig. 6a). Among Arthropoda, Maxillopoda were dominant at the class level; among them, Harpacti­ coida copepods were dominant (Fig. 6b). Among Annelida, all the OTUs were assigned to Polychaeta (Fig. 6d). The diversity indices calculated based on the number of reads for metazoan OTUs are shown in Table 6. Two-way PERMANOVA (factor “nucleotid”: DNA vs. RNA; factor “setting”: base vs. terrace) revealed that there were no significant differences in H0 or J’. On the contrary, there were significant differences between the “nucleotid” and “setting” factors in the number of OTUs, Expected, and Chao 1; however, there were no significant interactions. The plots based on the eDNA and eRNA data were plotted on the bottom left and top right sides of a PCO plot based on the present/absent data of metazoan OTUs (Fig. 7). Terrace and base stations were plotted on the top left and bottom right sides regardless of eDNA or eRNA. Twoway PERMANOVA (factor “nucleotid”: DNA vs. RNA; factor “setting”: base vs. terrace) revealed a significant difference in both factors, but there no interactions. SIMPROF analysis indicated that eDNA-based analysis could not distinguish differences in the settings, whereas eRNA-based analysis could.

3.2. Metabarcoding analysis MiSeq sequencing yielded a total of 1,588,250 reads, and 1,263,660 reads (average 361.1 bp in length) remained after quality check and chimera removal (Table S1). A total of 1,158,099 reads and 3886 OTUs were obtained after OTU construction, alignment, and singleton removal. Among those, 1,152,437 reads and 3862 OTUs were assigned to Eukaryota. The numbers of reads and OTUs assigned to eukaryotes and metazoans derived from the eDNA- and eRNA-based data are shown in Table 4. The number of reads from the eDNA and eRNA datasets were comparable, but there were more OTUs from eRNA than that from eDNA. The 3886 OTUs were classified as OTUs derived from the eDNA dataset only (305 OTUs), those from the eRNA dataset only (2064 OTUs), and those from the both the eDNA and eRNA datasets (shared OTUs, 1517 OTUs). Fig. 4 shows the relative taxonomic proportions obtained from these three datasets. Compared to the DNA-only dataset, the RNA-only dataset exhibited a higher proportion of fungi within Opisthokonta as well as Ciliophora within Alveolata. Venn diagrams showing the numbers of OTUs from the eDNA and eRNA datasets assigned to metazoans at each station are shown in Fig. 5. The propor­ tion of shared OTUs among all stations was approximately 30–80%. At stations other than B02, the proportions of RNA-only OTUs were greater than those of DNA-only OTUs.

4. Discussion 4.1. Evaluation of meiofaunal assemblage by conventional methods We failed to detect differences in meiofaunal density, diversity, and assemblage structures between stations in JA06 Seamount using the conventional methods. Wei et al. (2010) compiled meiofaunal densities worldwide and elucidated the relationship between meiofaunal density and water depth. According to their equation, the predicted meiofaunal densities at the base (4500 m) and terrace (1300 m) of JA06 Seamount are 147 and 369 ind./10 cm2, respectively. However, the meiofaunal densities observed in our study area were substantially lower than these estimated densities (Table 2). Shirayama (1984) reported comparable low meiofaunal densities in the abyssal plain of the Western Pacific Ocean as those measured in our study. This low meiofaunal density would be due to the low productivity at this study area, which caused difficulties in the detection of significant differences in meiofaunal as­ semblages. Nematoda and Copepoda are the most dominant taxa in the deep-sea environment, accounting for more than 90% (e.g. Giere, 2009),

Table 2 Meiofaunal density (ind./10 cm2) in the 0–5 cm sediment layer obtained from each station. B02 B02-2 B03 B04 T01 T02-2

Desmoscolecoidae

Other nematoda

Copepoda

Nauplii

Polychaeta

Bivalvia

Ostracoda

Kinorhyncha

Tardigrada

Gastrotricha

Total

0.0 1.9 0.0 0.0 2.1 0.4

2.6 13.0 2.8 16.1 44.4 61.5

0.2 3.5 0.7 3.7 10.7 5.5

1.2 4.9 0.2 3.0 9.8 3.5

0.0 0.0 0.0 0.2 0.7 0.7

0.0 0.0 0.0 0.2 0.0 0.0

0.2 0.5 0.0 1.2 0.9 0.8

0.0 0.2 0.0 0.0 0.2 0.1

0.0 0.5 0.0 0.0 0.0 0.0

0.0 0.0 0.0 0.0 0.0 0.1

4.2 24.4 3.7 24.4 68.9 72.5

4

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

Fig. 2. Vertical distribution of meiofauna at each station.

and the taxa except for these dominant taxa contribute to the difference in meiofaunal assemblages in the large-scale analysis (Bianchelli et al., 2010). Zeppilli et al. (2013) compared meiofaunal assemblages between the summit, flank, and base of Condor Seamount in the Northeast Atlantic Ocean by using the conventional microscopic methods and observed significant differences in meiofaunal density, biomass, and composition at high taxonomic levels at the depth of 200 m (seamount summit) to 1700 m (seamount base). The meiofaunal density at their study area exceeded 200 ind./10 cm2, which was higher than that observed at JA06 Seamount. In their study, the significant differences in the assemblages of meiofauna was detected due to the high meiofaunal density in the studied area. In this study, however, sufficient number of meiofaunal specimens could not be examined to detect significant dif­ ferences due to the low meiofaunal density at our study area.

4.2. Comparison of the eDNA- and eRNA-based approaches In the present study, eRNA-based analysis yielded more OTUs than eDNA-based analysis. This result was contrary to the expectation, because DNA was well preserved after cells die. Such “legacy DNA” in the sediments would be expected to lead to eDNA-based analysis yielding more OTUs. Guardiola et al. (2016) reported that eDNA-based analysis yielded more OTUs than eRNA-based analysis, and they also found that eRNA-based analysis yielded more OTUs assigned to Amoe­ bozoa and Nematoda than eDNA-based analysis. Most OTUs that were derived solely from eRNA-based analysis had few reads, indicating that they were not amplified by eDNA-based analysis (Guardiola et al., 2016). Here, in our study area, legacy DNA masked rare OTUs for all taxa, and eRNA-based analysis yielded more OTUs than eDNA-based analysis. In addition, at B02, eDNA-based analysis yielded more OTUs than eRNA-based analysis. Given that diatomaceous ooze was deposited at this station, living meiofauna would be scarce; therefore, legacy DNA 5

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

Fig. 3. PCO plot based on meiofaunal density at each station in (a) the 0–5 cm sediment layer and (b) the 0–1 cm sediment layer.

was amplified, thus masking the living OTUs at this station. eDNA- and eRNA-based analyses resulted in a total of 691 OTUs (275,522 reads) assigned to metazoan. Although these numbers cannot be simply compared to the meiofaunal density and number of taxa ob­ tained using the conventional methods, the metabarcoding analysis produced a large amount of data. Thus, the analysis succeeded in revealing significant differences in the spatial variation of meiofaunal diversity and community structure within the seamount, specifically between the terrace and base. Moreover, because eDNA and eRNA behave in a different manner in natural environments, the community structures determined by eDNA- and eRNA-based analyses were signif­ icantly different. Guardiola et al. (2016) observed a similar spatiotem­ poral pattern according to eDNA- and eRNA-based analyses, showing that both analyses yielded different community structures. Thus, legacy DNA did not obfuscate the patterns (Guardiola et al., 2016). Further­ more, Laroche et al. (2017, 2018) noted that eRNA-based analysis is less sensitive than eDNA-based analysis in detecting the effects of the oil and gas drilling activities. The amount of ribosomal RNA varies according to physiological conditions and even taxa, which would mask the actual pattern that can be detected by eDNA-based analysis (Laroche et al., 2017, 2018). In addition, Orsi et al. (2013) suggested that the ribosomal

Table 4 Numbers and proportions of the Reads and OTUs at each station. Eukaryota DNA B02 B02-2 B03 B04 T01 T02-2 RNA B02 B02-2 B03 B04 T01 T02-2 Total

Metazoa

Metazoa proportion, %

OTUs

Reads

OTUs

Reads

OTUs

Reads

1803 605 356 631 394 695 557 3524 790 890 1834 1642 2343 2256 3826

572,945 98,208 117,807 102,928 105,079 89,611 59,312 579,492 104,316 99,249 86,488 97,095 89,132 103,212 1,152,437

376 116 63 120 88 176 121 603 97 102 232 245 366 324 691

126,988 22,300 16,256 16,841 28,308 26,462 16,821 148,534 6332 7257 11,810 57,369 28,873 36,893 275,522

20.9 19.2 17.7 19.0 22.3 25.3 21.7 17.1 12.3 11.5 12.6 14.9 15.6 14.4 18.1

22.2 22.7 13.8 16.4 26.9 29.5 28.4 25.6 6.1 7.3 13.7 59.1 32.4 35.7 23.9

Fig. 4. Donut charts of OTUs for the DNA-only, shared eDNA/eRNA, and RNA-only datasets. The charts show the relative number of OTUs at the highest assigned taxonomic levels. The numbers of OTUs for each dataset are indicated inside each donut chart. 6

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

Fig. 5. Venn diagrams showing the number of the OTUs based on eDNA- and eRNA-based analyses at each station. Table 5 Numbers and proportions of the Reads and OTUs for each metazoan taxon. DNA Metazoa Porifera Mesozoa Cnidaria Annelida Arthropoda Bryozoa Chaetognatha Chordata Gastrotricha Hemichordata Mollusca Nematoda Nemertea Platyhelminthes Rotifera Kinorhyncha Loricifera Tardigrada Xenacoelomorpha uncultured

RNA

OTUs

Reads

376 2 0 186 23 57 1 0 3 13 1 1 55 1 12 1 0 2 0 15 3

126,988 186 0 59,121 12,035 26,970 135 0 162 2044 261 288 16,104 14 5882 99 0 42 0 3415 230

OTU % 0.5 0.0 49.5 6.1 15.2 0.3 0.0 0.8 3.5 0.3 0.3 14.6 0.3 3.2 0.3 0.0 0.5 0.0 4.0 0.8

RNA of some eukaryotic taxa might persist longer in the sediments than previously believed. Whether the eDNA- and eRNA based analyses show the same or different patterns of metazoan communities depends on the composition and activity of organisms at the study area. Nevertheless, if we assume that eDNA- and eRNA-based methods provide a reasonable picture of diversity despite several uncertainties, both methods are useful for such purposes because the true picture can be obtained only after a sophisticated morphological analysis is conducted.

Reads%

OTUs

Reads

OTU %

Reads%

0.1 0.0 46.6 9.5 21.2 0.1 0.0 0.1 1.6 0.2 0.2 12.7 0.0 4.6 0.1 0.0 0.0 0.0 2.7 0.2

603 2 1 283 26 73 2 1 3 20 0 2 130 1 24 0 2 5 1 24 3

148,534 17 96 35,496 49,349 23,304 46 31 300 4583 0 166 20,495 74 8572 0 75 1739 28 3355 808

0.3 0.2 46.9 4.3 12.1 0.3 0.2 0.5 3.3 0.0 0.3 21.6 0.2 4.0 0.0 0.3 0.8 0.2 4.0 0.5

0.0 0.1 23.9 33.2 15.7 0.0 0.0 0.2 3.1 0.0 0.1 13.8 0.0 5.8 0.0 0.1 1.2 0.0 2.3 0.5

focused on the OTUs (≒ species). However, the comparison would be possible at the higher taxonomic level. eDNA- and eRNA-based metabarcoding analyses can examine meiofaunal assemblages with only small sediment samples, and can provide large amounts of genetic information. In the present study, Cnidaria, Annelida, Arthropoda, and Nematoda were the dominant phyla in terms of numbers of reads and OTUs. Harpacticoid copepods accounted for a large percentage of arthropods. Microscopic observa­ tions showed that nematodes and copepods (mostly classified to Har­ pacticoida) were dominant in the deep-sea sediment, and this finding was corroborated by the literature (Giere, 2009). Furthermore, our microscopic observations were consistent with our metabarcoding re­ sults. All OTUs assigned to Annelida were assigned to Polychaeta, which is usually the third-most abundant taxon in deep-sea metazoan meio­ faunal assemblages after Nematoda and Harpacticoida (Giere, 2009). Accordingly, Polychaeta was the third-most abundant taxon in the

4.3. Comparison between the conventional and metabarcoding approaches The results obtained based on the conventional and metabarcoding methods could not be simply compared in this study because the two methods focused on the different taxonomic group: the former focused on the high taxonomic group (Nematoda, Copepoda, etc), the later 7

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx

Fig. 6. Donut charts of OTUs for the four dominant phyla: (a) Cnidaria, (b) Arthropoda, (c) Nematoda, and (d) Annelida. The numbers of OTUs assigned to each taxa are indicated inside each donut chart. Table 6 Metazoan diversity based on eDNA- and eRNA-based analyses at each station. Num of OTUs B02 B02-2 B03 B04 T01 T02-2

H0

J0

Expected

Chao 1

DNA

RNA

DNA

RNA

DNA

RNA

DNA

RNA

DNA

RNA

116 63 120 88 176 121

97 102 232 245 366 324

2.90 2.47 3.60 3.21 4.29 3.58

3.12 3.26 4.06 2.11 3.98 3.25

0.61 0.60 0.75 0.72 0.83 0.75

0.68 0.71 0.75 0.38 0.67 0.56

88.0 44.0 92.0 66.9 145.1 96.0

90.9 94.6 200.6 138.3 261.9 222.8

147.6 105.0 159.0 109.0 197.9 138.8

172.6 125.3 282.2 296.3 411.3 363.3

present study, although they exhibited considerably low density. Because polychaetes fragment easily during fixation, only complete specimens and fragments with heads were counted in the present study. Moreover, many polychaetes are referred to “temporary meiofauna” because they only fall within the meiofaunal size range during their juvenile stages. Juvenile specimens are difficult to be recognized and, therefore, cannot be counted under a microscope. Thus, the conven­ tional methods might have underestimated polychaete density. Approximately half of the total OTUs were assigned to Cnidaria both in the results obtained using eDNA- and eRNA-based metabarcoding (Table 5), although they were rare taxa in the deep-sea meiofaunal assemblage. Such a high abundance of OTUs assigned to Cnidaria has not been reported by metabarcoding studies in deep-sea environments (Bik et al., 2012; Guardiola et al., 2016; Sinniger et al., 2016). Among the Cnidaria OTUs, the most abundant class was Siphonophorae, which comprises typical planktonic organisms. The carcasses of such gelati­ nous zooplankton have been reported to be deposited on the deep-sea floor (Billett et al., 2006). Cnidaria OTUs were detected not only by eDNA-based analysis but also by eRNA-based analysis (Table 5), which did not mean that the legacy DNA of Cnidaria remained without

undergoing decomposition. These results suggested that the organic matter of carcasses of these planktonic organisms deposited in the deep-sea zone below 1000 m depth, including the 4500 m-deep JA06 Seamount, decomposed quite slowly. It is worth mentioning that “soft-bodied” meiofauna (e.g., Gastro­ tricha, Platyhelminthes, Rotifera, and Xenacoelomorpha) were detected by the eDNA- and eRNA-based analyses, but not by the conventional microscopic analysis (Table 5). These taxa readily break down after formalin fixation (Pfannkuche and Thiel, 1988), and knowledge on their taxonomy, diversity, and distribution is severely limited (Curini-Galletti et al., 2012). Thus, the present study showed that molecular and met­ abarcoding approaches was able to resolve the hidden diversity of these soft-bodied meiofauna, as suggested in previous studies (Leasi and Norenburg, 2014; Rzeznik-Orignac et al., 2017). 4.4. Concluding remarks Although seamounts are recognized as important areas in deep-sea ecosystems, there are growing concerns regarding the effects of anthropogenic activities, such as deep-sea trawling and resource 8

T. Kitahashi et al.

Deep-Sea Research Part I xxx (xxxx) xxx illuminates biogeographic patterns in deep-sea and shallow water sediments. Mol. Ecol. 21, 1048–1059. Billett, D.S.M., Bett, B.J., Jacobs, C.L., Rouse, I.P., Wigham, B.D., 2006. Mass deposition of jellyfish in the deep Arabian Sea. Limnol. Oceanogr. 51, 2077–2083. Caporaso, J.G., Bittinger, K., Bushman, F.D., DeSantis, T.Z., Andersen, G.L., Knight, R., 2010a. PyNAST: a flexible tool for aligning sequences to a template alignment. Bioinformatics 26, 266–267. Caporaso, J.G., Kuczynski, J., Stombaugh, J., Bittinger, K., Bushman, F.D., Costello, E.K., Fierer, N., Pe~ na, A.G., Goodrich, J.K., Gordon, J.I., Huttley, G.A., Kelley, S.T., Knights, D., Koenig, J.E., Ley, R.E., Lozupone, C.A., McDonald, D., Muegge, B.D., Pirrung, M., Reeder, J., Sevinsky, J.R., Turnbaugh, P.J., Walters, W.A., Widmann, J., Yatsunenko, T., Zaneveld, J., Knight, R., 2010b. QIIME allows analysis of highthroughput community sequencing data. Nat. Methods 7, 335–336. Chao, A., Chazdon, R.L., Colwell, K.R., Shen, T.J., 2005. A new statistical approach for assessing similarity of species composition with incidence and abundance data. Ecol. Lett. 8, 148–159. Clark, M.R., Rowden, A.A., Schlacher, T., Williams, A., Consalvey, M., Stocks, K.I., Rogers, A.D., O’Hara, T.D., White, M., Shank, T.M., Hall-spencer, J.M., 2010. The ecology of seamounts: structure, function, and human impacts. Annu. Rev. Mar. Sci. 2, 253–278. Corinaldesi, C., Barucca, M., Luna, G.M., Dell’Anno, A., 2011. Preservation, origin and genetic imprint of extracellular DNA in permanently anoxic deep-sea sediments. Mol. Ecol. 20, 642–654. Cristescu, M.E., 2019. Can environmental RNA revolutionize biodiversity science? Trend. Evol. Ecol. 34, 694–697. Curini-Galletti, M., Artois, T., Delogu, V., de Smet, W.H., Fontaneto, D., Jondelius, U., Leasi, F., Martínez, A., Meyer-Wachsmuth, I., Nilsson, K.S., Tongiorgi, P., Worsaae, K., Todaro, M.A., 2012. Patterns of diversity in soft-bodied meiofauna: dispersal ability and body size matter. PloS One, e33801. Danovaro, R., 2010. Methods for the Study of Deep-Sea Sediments, Their Functioning and Biodiversity. CRC Press. Dell’Anno, A., Danovaro, R., 2005. Extracellular DNA plays a key role in deep-sea ecosystem functioning. Science 309, 2179. George, K.H., 2013. Faunistic research on metazoan meiofauna from seamounts-a review. Meiofauna Marina 20, 1–32. George, K.H., Pointner, K., Packmor, J., 2018. The benthic Copepoda (Crustacea) of anaximenes seamount (eastern mediterranean sea)—community structure and species distribution. Prog. Oceanogr. 165, 299–316. Giere, O., 2009. Meiobenthology, second ed. Springer, Berlin. Guardiola, M., Wangensteen, O.S., Taberlet, P., Coissac, E., Uriz, M.J., Turon, X., 2016. Spatio-temporal monitoring of deep-sea communities using metabarcoding of sediment DNA and RNA. PeerJ 4, e2807. Hurlbert, H.S., 1971. The nonconcept of species diversity: a critique and alternative parameter. Ecology 52, 577–586. Laroche, O., Wood, S.A., Tremblay, L.A., Ellis, J.I., Lear, G., Pochon, X., 2018. A crosstaxa study using environmental DNA/RNA metabarcoding to measure biological impacts of offshore oil and gas drilling and production operations. Mar. Pollut. Bull. 127, 97–107. Laroche, O., Wood, S.A., Tremblay, L.A., Lear, G., Ellis, J.I., Pochon, X., 2017. Metabarcoding monitoring analysis: the pros and cons of using co-extracted environmental DNA and RNA data to assess offshore oil production impacts on benthic communities. PeerJ 5, e3347. Leasi, F., Norenburg, J.L., 2014. The necessity of DNA taxonomy to reveal cryptic diversity and spatial distribution of meiofauna, with a focus on nemertea. PloS One 9, e104385. Lejzerowicz, F., Voltsky, I., Pawlowski, J., 2013. Identifying active foraminifera in the Sea of Japan using metatranscriptomic approach. Deep-Sea Res. II 86–87, 214–220. Markmann, M., Tautz, D., 2005. Reverse taxonomy: an approach towards determining the diversity of meiobenthic organisms based on ribosomal RNA signature sequences. Philos. Trans. R. Soc. B 360, 1917–1924. Martin, M., 2011. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet Journal 17, 10–12. Oksanen, J., Guillaume Blanchet, F., Friendly, M., Kindt, R., Legendre, P., McGlinn, D., Minchin, P.R., O’Hara, R.B., Simpson, G.L., Solymos, P., Henry, M., Stevens, H., Szoecs, E., Wagner, H., 2018. Vegan: community ecology package. R package version 2.5-3. https://CRAN.R-project.org/package¼vegan. Orsi, W., Biddle, J.F., Edgcomb, V., 2013. Deep sequencing of subseafloor eukaryotic rRNA reveals active fungi across marine subsurface provinces. PloS One 8, e56335. Pawlowski, J., Christen, R., Lecroq, B., Bachar, D., Shahbazkia, H.R., Amaral-Zettler, L., Guillou, L., 2011. Eukaryotic richness in the abyss: insights from pyrotag sequencing. PloS One 6 (4), e18169. Pfannkuche, O., Thiel, H., 1988. Sample processing. In. In: Higgins, R.P., Thiel, H. (Eds.), Introduction to the Study of Meiofauna. Smithsonian Institute Press, Washington, D. C., pp. 134–145 Pochon, X., Zaiko, A., Fletcher, L.M., Laroche, O., Wood, S.A., 2017. Wanted dead or alive? Using metabarcoding of environmental DNA and RNA to distinguish living assemblages for biosecurity applications. PloS One 12, e0187636. Probert, P.K., Christiansen, S., Gjerde, K.M., Gubbay, S., Santos, R.S., 2007. Management and conservation of seamounts. In: Pitcher, T.J., Morato, T., Hart, P.J.B., Clark, M.R., Haggan, N., Santos, R.S. (Eds.), Seamounts: Ecology, Fisheries and Conservation. Blackwell Publishing, Oxford, UK, pp. 442–475. Quast, C., Pruesse, E., Yilmaz, P., Gerken, J., Schweer, T., Yarza, P., Peplies, J., Gl€ ockner, F.O., 2013. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 41, 590–596. Rowden, A.A., Schlacher, T.A., Williams, A., Clark, M.R., Stewart, R., Althaus, F., Bowden, D.A., Consalvey, M., Robinson, W., Dowdney, J., 2010. A test of the

Fig. 7. PCO plot based on the present/absent data in the OTUs according to the eDNA- and eRNA-based analyses. Circles indicate significant groups (SIMPROF, p < 0.05).

exploitation, on these ecosystems (Probert et al., 2007; Clark et al., 2010). The present study suggested that metabarcoding analysis using eDNA or eRNA was able to provide invaluable information on meio­ faunal assemblages in seamounts even when the conventional micro­ scopic analysis provided insufficient information to elucidate the heterogeneity of meiofaunal assemblages within a seamount due to low meiofaunal density. Declaration of competing interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgements The authors are grateful to the officers, crew, and the Chief Scientist, Dr. Masashi Hayakawa of the R/V Hakurei cruise. Special thanks go to Dr. Kazuya Naito, Dr. Hideki Sugishima (JOGMEC), Dr. Hiroshi Ama­ kawa, Dr. Yuya Tada (JAMSTEC), and Dr. Yohei Taketomo (Ocean En­ gineering & Development Co., Ltd.) who helped with sampling onboard. We would like to thank to the two anonymous reviewers who provided helpful comments on an earlier draft of the manuscript. Generic Map­ ping Tools (GMT, Wessel and Smith, 1995) was used in this study. This research was partly supported by Cross-ministerial Strategic Innovation Promotion Program (SIP) from Cabinet Office (CAO), Government of Japan. This research was conducted as part of the commissioned projects by Agency for Natural Resources and Energy, Ministry of Economy, Trade and Industry. Appendix A. Supplementary data Supplementary data to this article can be found online at https://doi. org/10.1016/j.dsr.2020.103253. References Barnett, P.R.O., Watson, J., Connelly, D., 1984. A multiple corer for taking virtually undisturbed samples from shelf, bathyal and abyssal sediments. Oceanol. Acta 7, 399–408. Bianchelli, S., Gambi, C., Zeppilli, D., Danovaro, R., 2010. Metazoan meiofauna in deepsea canyons and adjacent open slopes: a large-scale comparison with focus on the rare taxa. Deep-Sea Res. II 57, 420–433. Bik, H.M., Sung, W., De Ley, P., Baldwin, J.G., Sharma, J., Rocha-Olivares, A., Thomas, W.K., 2012. Metagenetic community analysis of microbial eukaryotes

9

Deep-Sea Research Part I xxx (xxxx) xxx

T. Kitahashi et al.

Wei, C.-L., Rowe, G.T., Escobar-Briones, E., Boetius, A., Soltwedel, T., Caley, M.J., Soliman, Y., Huettmann, F., Qu, F., Yu, Z., Pitcher, C.R., Haedrich, R.L., Wicksten, M. K., Rex, M.A., Baguley, J.G., Sharma, J., Danovaro, R., MacDonald, I.R., Nunnally, C. C., Deming, J.W., Montagna, P., L�evesque, M., Weslawski, J.M., WlodarskaKowalczuk, M., Ingole, B.S., Bett, B.J., Billett, D.S.M., Yool, A., Bluhm, B.A., Iken, K., Narayanaswamy, B.E., 2010. Global patterns and predictions of seafloor biomass using random forests. PloS One 5, e15323. Wessel, P., Smith, W.H., 1995. New version of the generic mapping tools released. EOS Trans. AGU 76, 329. Yamamoto, H., Hashimoto, Y., Ezaki, T., 1996. Study of nonculturable Legionella pneumophila cells during multiple-nutrient starvation. FEMS Microbiol. Ecol. 20, 149–154. Zeppilli, D., Bongiorni, L., Cattaneo, A., Danovaro, R., Santos, R.S., 2013. Meiofauna assemblages of the condor seamount (North-east Atlantic ocean) and adjacent deepsea sediments. Deep-Sea Res. II 98, 87–100. Zhang, J., Kobert, K., Flouri, T., Stamatakis, A., 2014. PEAR: a fast and accurate Illumina Paired-End reAd mergeR. Bioinformatics 30, 614–620.

seamount oasis hypothesis: seamounts support higher epibenthic megafaunal biomass than adjacent slopes. Mar. Ecol. 31, 95–106. RStudio Team, 2016. RStudio. Integrated Development for R. RStudio, Inc., Boston, MA. http://www.rstudio.com/. Rzeznik-Orignac, J., Kalenitchenko, D., Mariette, J., Bodiou, J.-Y., Le Bris, N., Derelle, E., 2017. Comparison of meiofaunal diversity by combined morphological and molecular approaches in a shallow Mediterranean sediment. Mar. Biol. 164, 40. Sautya, S., Ingole, B., Ray, D., St€ ohr, S., Samudrala, K., Raju, K.K.A., Mudholkar, A., 2011. Megafaunal community structure of Andaman seamounts including the BackArc Basin - a quantitative exploration from the Indian Ocean. PloS One 6, e16162. Shirayama, Y., 1984. The abundance of deep-sea meiobenthos in the western pacific in relation to environmental factors. Oceanol. Acta 7, 113–121. Sinniger, F., Pawlowski, J., Harii, S., Gooday, A.J., Yamamoto, H., Chevaldonn�e, P., Cedhagen, T., Carvalho, G., Creer, S., 2016. Worldwide analysis of sedimentary DNA reveals major gaps in taxonomic knowledge of deep-sea benthos. Front. Mar. Sci. 3, 92. Takishita, K., Kakizoe, N., Yoshida, T., Maruyama, T., 2010. Molecular evidence that phylogenetically diverged ciliates are active in microbial mats of deep-sea cold-seep sediment. J. Eukaryot. Microbiol. 57, 76–86.

10