Skip to main content

Comparison of the transcriptomes of two tardigrades with different hatching coordination

Abstract

Background

Tardigrades are microscopic organisms, famous for their tolerance against extreme environments. The establishment of rearing systems of multiple species has allowed for comparison of tardigrade physiology, in particular in embryogenesis. Interestingly, in-lab cultures of limnic species showed smaller variation in hatching timing than terrestrial species, suggesting a hatching regulation mechanism acquired by adaptation to their habitat.

Results

To this end, we screened for coordinated gene expression during the development of two species of tardigrades, Hypsibius exemplaris and Ramazzottius varieornatus, and observed induction of the arthropod molting pathway. Exposure of ecdysteroids and juvenile hormone analog affected egg hatching but not embryonic development in only the limnic H. exemplaris.

Conclusion

These observations suggest a hatching regulation mechanism by the molting pathway in H. exemplaris.

Background

Tardigrades are microscopic animals found in various meiofauna and are prominently known for their capability to tolerate various extreme environments during the ametabolic state designated as cryptobiosis [23]. Recent advances in the establishment of tardigrade genomic information has allowed for comprehensive molecular analysis in tardigrade species [17, 41, 47, 48].

The establishment of rearing systems for several species has enabled various quantitative zoological and molecular analyses, mainly focused on the limnic Hypsibius exemplaris and terrestrial Ramazzottius varieornatus [14, 18]. In particular, developmental research of the H. exemplaris embryo has provided valuable evidence regarding the evolution and phylogeny of tardigrades within Ecdysozoa (Reviewed in Altiero et al. [3]). A large amount of developmental research has focused on the morphology of H. exemplaris, due to their rapid embryogenesis and life cycle. The H. exemplaris embryo hatches after around 4–4.5 days [14]. Observations based on microscopy-based techniques have reported that the first cleavage is observed at 2 h post laying (hpl), and gastrulation is completed by 16.5 hpl [14]. Further Scanning Electron Microscope (SEM) observations indicate that the head and limb buds, segmental lines, claws and mouth are visible by − 25 hpl, − 30 hpl and 35 and 50 hpl, respectively [16], suggesting that embryonic development is concluded around 3 days after oviposition. On the other hand, establishment of gene sequence databases and genetic engineering methods has enabled the molecular analysis of development related genes [9, 17, 48]. Immunostaining-based expression analysis and in situ hybridization analysis has enabled the analysis of expression patterns of several genes: pair-rule gene paired (Pax3/7), the segment polarity gene engrailed [13], and Homeobox genes [39]. Additionally, RNA interference knock-down of the expression of mago-nashi resulted in failure of elongation in the developing embryo, suggesting a major role of mago-nashi during embryogenesis [43].

On the other hand, the hatching timing of in-lab cultures of R. varieornatus extended to 5.8 ± 1.2 days (mean ± SD), slightly longer than H. exemplaris [18]. Relatively higher diversity in hatching timing of in-lab cultures which lack non-predictable stimuli (i.e. environmental factors) suggests that the diversity may be the result of adaptation to their living habitat.

To this end, we sought to identify factors that cause diversity in hatching timing. We have previously analyzed the genomes of two tardigrades from different habitats, the limnic H. exemplaris and terrestrial R. varieornatus. Genomes of both species have been sequenced, making these two species ideal for a comparative genomic/transcriptomic analysis, such as the identification of anhydrobiosis related genes [17, 48]. We first observed hatching timing in H. exemplaris to calculate the diversity in this species. We then conducted a comparative transcriptomic analysis of the embryo and juvenile stages in both species to identify pathways that are induced around 1 day prior to hatching. We finally exposed embryos to the ecdysteroid 20-Hydroxyecdysone (20-E) and juvenile hormone analog Fenoxycarb to understand the effects of such substances on the embryo. We observed that H. exemplaris also has tightly regulated embryo hatching, and the arthropod molting pathway was induced. 20-E and Fenoxycarb exposure induced a latent state in only the H. exemplaris embryo. These findings support the idea that limnic species may have tighter regulation in embryo hatching and that embryo hatching may be hormone regulated.

Results

Previous studies have shown that embryonic development occurs in 4–4.5 days in H. exemplaris, and generation time takes 13–14 days. We collected embryos immediately after oviposition, and observed the development time, as well as generation time (Fig. 1a). Embryonic development took 4.03 ± 0.50 days from oviposition to hatch. The 1st and 2nd oviposition occurred around 7–11 days and 16–22 days after hatching, with 2.68 ± 0.75 and 8.34 ± 4.12 eggs per individual, respectively (Fig. 1b). The body length reached the average length of an adult (240 μm) at around 15d after oviposition (Additional file 1: Figure S1). In comparison, R. varieornatus requires 5.72 ± 1.13 days to hatch and the 1st and 2nd oviposition occurs at 8–12 days and 13–18 days after hatching. Additionally, we observed that the average number of eggs per individual was 7.37 ± 2.83, similar to 7.8 in previous studies [18]. We therefore assume that hatching related genes may be upregulated around Egg 3d and 5d in H. exemplaris and R. varieornatus, respectively.

Fig. 1
figure1

Comparison of hatching timing between H. exemplaris and R. varieornatus. a Probability density plot of the days required for egg hatching in both species. Seventy-three and 259 embryos were observed for H. exemplaris and R. varieornatus, respectively. A normal distribution curve is overlaid for each species. b Number of eggs that were oviposit in each day per individual in H. exemplaris and R. varieornatus. Twenty-one and 40 individuals were observed for each species, respectively. Data from Horikawa et al. [18] was used for R. varieornatus

We have previously sequenced the transcriptome of the embryo and juvenile stages in H. exemplaris, and, in part, R. varieornatus. Hence, we additionally sequenced the juvenile stages of R. varieornatus, and conducted a comprehensive analysis of developmental and juvenile stages of both species (Additional file 9: Table S1). High correlation was observed between replicates (N = 3). In order to validate our sequencing data, we first compared the expression profiles of HOX genes observed previously by single embryo RNA-Seq methods [26]. In our data, HOX genes showed an Egg 2-3d specific expression in H. exemplaris, while a wide expression was observed at the Egg 2-4d in R. varieornatus (Additional file 2: Figure S2). In particular, caudal was expressed in Egg 2d, while other HOX genes were expressed at Egg 2-3d. In comparison, caudal was mainly expressed around 500 min (8.3 h), while other HOX genes (i.e. lab, ftz, disformed, zen) were expressed around 1600 min (26.7 h) in the CEL-Seq data (Additional file 3: Figure S3). To obtain a general view of the transcriptome profile of both species, we conducted SOM clustering of gene expression profiles. We found two out of six clusters showed an Egg 3d-specific expression in H. exemplaris (Fig. 2, groups 1 and 2, Additional file 10: Data S1), suggesting dynamic shift in the embryonic transcriptome prior to hatching in only H. exemplaris. No such shift was observed in R. varieornatus.

Fig. 2
figure2

SOM clustering of gene expression profiles during development. SOM clustering of TPM values were performed in R. Genes with TPM values over 1 were used. The areas colored in blue and red represent embryo and juvenile stages, respectively. The Y axis indicates Z-scale normalized TPM values. The arrows indicate an increase in Egg 3d of H. exemplaris

In order to identity potential mediators of hatching, we screened the 5436 genes induced at the Egg 3d in H. exemplaris (Group 1 and 2, Additional file 11: Data S2). It is commonly known that many factors are regulated during embryonic development; however, most of these mechanisms have yet to be validated in tardigrades. We focused on signaling through neuropeptides and hormones, of which we identified 232 receptor genes in H. exemplaris with gene ontology annotations containing “hormone”. Eighty-two of these orthologs were contained in groups 1 and 2. We obtained 20 orthologs that were significantly expressed in only Egg 3d. Although a majority of these genes were atrial natriuretic peptide receptors, we identified multiple orthologs of the arthropod molting pathway (i.e. EcR, RXR, E75, HR3, HR4). Therefore, we collected other genes in the molting pathway identified in previous studies and observed their expression patterns during development. The expression patterns of these genes were highly upregulated in Egg 3d in H. exemplaris and constitutively moderately between Egg 2-5d in R. varieornatus (Fig. 3).

Fig. 3
figure3

Gene expression patterns for molting pathways during development. Genes identified in Schumann et al. [37] were identified by BLASTP searches against amino acid sequences of H. exemplaris and R. varieornatus, and Z-scaled TPM values were visualized. E: Egg, B: Juvenile, act/tun: Adult stages. The identified orthologs are as follows; sad (OQV21685.1, g2400.t1), EcR (OQV14677.1, OQV23446.1, g12979.t1, g8354.t1), USP/RXR (OQV18794.1, OQV18795.1, g11953.t1, g11953.t2), E74 (OQV12233.1, g10224.t1), E75 (OQV18927.1, g7221.t1), E78 (OQV14742.1, g11712.t1), HR3 (OWA50673.1, g6536.t1), HR4 (OQV14243.1, g2273.t1), bFtz-F1 (OQV18443.1, g4792.t1), CYP18A1 (OQV16526.1, g13718.t1). An additional copy of EcR was found during this analysis (H. exemplaris OQV23446.1, R. varieornatus g8354.t1)

We then validated whether disrupting the molting pathway inhibits hatching in H. exemplaris, we exposed H. exemplaris embryos to chemicals (20-E: 8.42 μM and 100 μM, Fenoxycarb: 3.3 μM and 100 μM) and observed whether they affected embryo development. While we observed that approximately 90% of the embryos exposed to low concentrations hatched within 5 days of oviposition (Additional file 6: Figure S6), there were several cases where all of the embryos in the same clutch did not hatch (Additional file 7: Figure S7). We observed claw formation at this stage, which is consistent with previous body developmental observations. Washing these embryos to remove chemical exposure induced hatching in 57.1% (4/7, five E0 embryos and two E3 embryos, of which three E0 embryos died) and 71.4% (5/7, one F0 and F1 embryos, three F2 embryo and two F3 embryos, of which one F0 and F2 embryos died) of the embryos exposed to 20-E and Fenoxycarb, respectively. Exposure of R. varieornatus embryo to low concentrations did not show any effect in hatching (Additional file 7: Figure S7). On the contrary, exposure to high concentrations of Fenoxycarb significantly suppressed hatching in embryos exposed at Egg 0-2d (Fig. 4, Additional file 8: Figure S8, ANOVA and Tukey HSD, FDR < 0.05). For the embryos that did not hatch, we also observed claw formation and movement inside the eggshells regardless of exposure time, suggesting that the embryos had finished development and did not die (Additional file 13: Data S4). Exposure to the EcR antagonist cucurbitacin B had no effect on embryo hatching (Data not shown).

Fig. 4
figure4

High concentration fenoxycarb exposure inhibited hatching in H. exemplaris. Percentages of successful hatching embryos exposed in 100 μM Fenoxycarb or 100 μM 20-E until hatch. Embryo exposed to 0.03% ethanol in Volvic mineral water at Egg 0d were used as controls. Hatching ratio for each clutch are overlaid on the boxplot. The x-axis indicates the day to start exposing the chemicals after oviposition. High concentration treatment of Fenoxycarb from 0~2 days after oviposition significantly impaired hatching, whereas treatment at 3 days after oviposition and treatment with 20-E did not. ANOVA and Tukey HSD, * = FDR < 0.05, ** = FDR < 0.01, *** = FDR < 0.001

Additional file 13. Movie of embryos exposed to long-term high concentration of chemicals. The shells of embryo exposed to more than 24 h to high concentration of fenoxycarb were peeled off, and movement were observed. https://github.com/abs-yy/Yoshida_and_Sugiura_etal_2018/blob/master/ Sdata_7_highConc-long-movement.mp4.

Discussion

Previous studies on the life history of H. exemplaris have observed that the embryo hatches around 4–4.5 days after oviposition [14]. We have observed that the H. exemplaris embryo hatches in a strictly controlled manner within a short time range (4.03 ± 0.50 days) compared to the wider range observed in R. varieornatus (5.72 ± 1.13 day, see [18]). The small variance in H. exemplaris is similar to the limnic species Isohypsibius myrops [19]. It is worth noting that the R. varieornatus strain YOKOZUNA-1 was derived from a single individual [18], suggesting that genetic variation would not be higher than H. exemplaris, and therefore the hatching time variability is not due to the genetic variability. Therefore, we hypothesized that transcriptome regulation affects hatching timing, and factors that are related to hatching regulation may be induced around a day prior to hatching; day 3 in H. exemplaris and day 5 in R. varieornatus.

To understand the molecular pathways regulated prior to hatching, we analyzed our previously sequenced transcriptome data of the developmental stages of two tardigrade species, H. exemplaris and R. varieornatus, with newly generated data for the juvenile stages of R. varieornatus [9, 48]. First, we validated our sequencing data by comparison of the expression profiles of developmental genes (i.e. Homeobox genes) assayed in previous single embryo sequencing data and in-situ hybridization experimental data [26, 39]. In our data, caudal and engrailed expressed earlier than HOX orthologs, at around 24–48 h post laying, and other orthologs were expressed around 48-72 h, approximately 36 h prior to hatching. These data are similar to those of previous reports of HOX expression profiles, supporting the validity of our data. In comparison, R. varieornatus HOX orthologs were expressed around Egg 3d-4d, 24 h later than H. exemplaris. This period is also approximately 36 h prior to hatching [18]. Clustering based analysis of the embryo-juvenile time course indicated an abrupt shift in transcriptome profile at Egg 3d in H. exemplaris. The Egg 3d stage is immediately before the hatch in H. exemplaris (4.03 ± 0.53 days), which suggests that H. exemplaris may have some sort of trigger that drastically changes the transcriptome profile for embryo hatching. The lack of such a shift in R. varieornatus remains a surprise, as we have sequencing data until approximately 20 h prior to hatching. This may be due to the low quality of the sequencing data, observed by the low mapping ratio, or pathways induced in H. exemplaris may be expressed less abruptly, but rather gradually during the last few days of embryogenesis.

Analysis of genes induced at Egg 3d of H. exemplaris revealed the arthropod molting pathway was upregulated at Egg 3d. Other pathways are stated in Additional file 14 (with accompanying Additional file 4: Figure S4, Additional file 5: Figure S5 and Additional file 12: Data S3). Molting is an essential event for the development of animals in the Ecdysozoa, however, the fact that the arthropod molting pathway is induced in embryos of Daphnia magna, Bombyx mori and Drosophila melanogaster [11, 21, 38, 40], suggests a mechanism conserved within Arthropoda that may be related to embryo hatching. In ecdysozoans that have been investigated, the molting pathway is initiated by the signaling hormone Ecdysone, synthesized by the Halloween genes in hexapods. It is released as the active state 20-hydroxyecdysone (20-E, [20]). The released hormone binds to the ecdysone receptor (EcR) and ultraspiracle (USP) heterodimer, which regulates the expression of the early genes (HR3, HR4, E75, Ftz-F1). These respond and transmits the signal to the late response genes, leading to new cuticle biosynthesis. Although most of the downstream genes of EcR/USP were conserved, only shadow (sad) gene in the 20-E biosynthesis pathway were found in both H. exemplaris and R. varieornatus [37]. These genes were induced at Egg 3 in H. exemplaris and relatively constitutively expressed (most expressed at Egg 5d) in R. varieornatus, approximately 1 day prior to egg hatching. We did not observe a lag between early and late genes since we sampled the embryos every 24 h. Interestingly, these genes were not expressed prior to the initial molting around Juvenile 7d, whereas they were in R. varieornatus (Juvenile 4d). Lack of expression during molting in H. exemplaris may imply that this pathway may not participate in the H. exemplaris molting, but in other pathways. On the other hand, the effects of metamorphosis negative regulator juvenile hormones (JH) are well studied [22]. The downstream gene of JH signaling Krüppel homolog-1 (Kr-h1) negatively regulates E75 (Li et al., 2018). Additionally, expression of Kr-h1 inhibits ecdysone biosynthesis in prothoracic glands, where 20-E precursor ecdysone are produced. These suggests that juvenile hormones inhibit the arthropod molting pathway in both ecdysone production and signal transition. Fenoxycarb, a versatile substance that functions non-dependently against juvenile hormones that are species specific, has been used intensively as a JH analog for pathway inhibition [15, 31]. Although the conservation of signaling pathway of JH have not been identified for tardigrades, fenoxycarb functions as a JH analog for a variety of metazoan species, which implies that this may also have an effect on tardigrades. Similarly, we hypothesized that 20-E may function in tardigrades considering that the EcR/USP genes were conserved. To validate if the molting pathway regulates hatching behavior, we evaluated the effect of exposure to the ecdysteroid 20-E and juvenile hormone analog Fenoxycarb. Exposure to high concentrations of Fenoxycarb or 20-E induced un-hatched embryos, suggesting that hormone exposure may have affected embryonic development or hatching. We confirmed that the majority of the un-hatched embryos hatched after hormone removal developed claws, suggesting that exposure to 20-E and fenoxycarb regulates only hatching and not embryonic development. On the contrary, exposure of both chemicals to R. varieornatus did not affect hatching timing. This was not the result of the lack of expression of related genes; the expression of EcR/USP downstream genes were constitutively expressed at moderate levels in the R. varieornatus embryo (higher levels prior to the initial molting). These suggests that similar hatching regulation mechanisms may not be conserved in R. varieornatus. Additionally, previous studies have suggested that tardigrade embryos may enter a “rest state” induced by environmental perturbations [1, 2]. Although there have been no reports of a rest state induced in H. exemplaris, hatching of hatch-inhibited embryo after the removal of hormones is remarkably similar. These observations suggest that the embryo hatching process can be regulated by chemicals, both positively and negatively. The differences between the time of effects between chemicals may be due to the differences of incorporation rates or the expression of receptors. Observation of inhibition of embryo hatch in high concentration H. exemplaris exposure suggests that a portion of the chemicals are penetrating the eggshell, of which such concentrations do not affect development but only the hatch regulation mechanism.

Together, these findings suggest that the variation in embryo hatching may suggest an adaptation against their living habitat; xerophilic R. varieornatus inhabits a desiccation-prone environment, and hygrophilic H. exemplaris is less likely to face environmental desiccation on a regular basis. Therefore, R. varieornatus may possess variation in hatching timing in order to avoid loss of entire broods that faces environmental stress close to hatching, therefore hatching regulation by the molting pathway may not be required. On the contrary, H. exemplaris may allow for controlled hatching regulation through internal hormones to enhance proliferation by a more rapid development. However, we have yet to validate if the phenotypic changes to the embryo originates from the chemical exposure; the lack of phenotypic change in low concentration exposure in both species may be result of failure of chemicals to penetrate the eggshell or 20-E/Fenoxycarb may not be the proper chemicals to use as juvenile hormone analogs. Furthermore, our findings are based on only two species from two families, Hypsibiidae and Ramazzottiidae; molecular data on other species or families will be required to validate our hypothesis. We await further studies to confirm if these findings are general to the Tardigrada phylum.

Conclusions

We conducted a comparative transcriptomic study using two tardigrade species to identify pathways induced for hatching regulation, and then validated whether these critical points could be interfered with by chemical exposure. We observed that the 3rd day of the embryo in H. exemplaris has a unique transcriptome profile compared to other stages, approximately 1 day prior to hatching, presumably for the strict control of its timing. Mirroring this observation, exposure of eggs to chemicals before the 3rd day influenced hatching only in H. exemplaris, which suggests that related pathways may be a critical component of H. exemplaris embryo hatching.

Methods

Tardigrade rearing and life history analysis

R. varieornatus was maintained using the methods previously described [18]. Briefly, individuals were placed on 2% Bacto agar (Difco) gels prepared with Volvic mineral water and fed with Chlorella vulgaris (Chlorella Industry). Individuals were moved to a new plate every week. The cultures were examined daily to obtain eggs and were observed every day to obtain juvenile samples for days 3–7.

H. exemplaris was also maintained using the same method with minor modifications (1.2% Bacto agar gel plates). For the life history analysis, 21 hatchlings were isolated from the culture plates, and were imaged every day. Body length was measured using ImageJ, and the days required for vitellogenesis, laying eggs days, and the number of eggs in each clutch were noted. Seventy-three eggs were used for hatching day observation. Standard variation was calculated for the hatching timing.

Transcriptome sequencing and data preparation

Isolated R. varieornatus specimens were collected and were subjected to transcriptome sequencing using the methods previously described with three replicates [9]. In brief, specimens were thoroughly washed with Milli-Q water on a sterile nylon mesh (Millipore), and were placed in a low-binding PCR tube. Pipette tips were used to homogenize the animals and lysed with the TRIzol reagent (Life Technologies). Total RNA was extracted using a Direct-zol RNA kit (Zymo Research), and the quality was validated with High Sensitivity RNA ScreenTape on TapeStation 2200 (Aglient Technologies). The mRNA was amplified with the SMARTer Ultra Low Input RNA Kit for Sequencing v.3, and Illumina libraries were prepared with the KAPA HyperPlus Kit (KAPA Biosystems). The library was quantified using the Qubit Fluorometer (Life Technologies), and size distribution was validated with TapeStation D1000 ScreenTape (Agilent Technologies). The library was size selected for 200 bp by manually cutting out an agarose gel and was purified with a NucleoSpin Gel and PCR Clean-up kit (Clontech). The samples were sequenced using the NextSeq 500 High Output Mode 75 cycles kit (Illumina) as single end reads. Adaptor sequences were removed, and the reads were de-multiplexed using the bcl2fastq v. 2 software (Illumina). All mRNA-Seq data have been uploaded to GEO under the accession ID GSE130111.

Informatics analysis

We previously performed time-series transcriptome analysis using the developmental stages of both H. exemplaris and R. varieornatus (Egg 1-5d), in addition to the juvenile stages of H. exemplaris (Juvenile 1-7d) and R. varieornatus (Juvenile 1d) [48]. These RNA-Seq data contained 4~17 M reads per sample. mRNA-Seq and small RNA-Seq data were downloaded from SRP098585 and SRX2495676, respectively. Sequence data were validated for sequencing quality using the FastQC software [5].

Amino acid and coding sequences and annotations for R. varieornatus and H. exemplaris were downloaded from www.ensembl.tardigrades.org [48]. Gene expression values were quantified using Kallisto v. 0.42.4 with the options: --bias -b 100 -t 31 --single -l 350 -s 50 [12]. The expression profiles were clustered by Self-Organizing Map (SOM, 2 × 3) using transcripts with an average over 1 TPM in all samples. To identify differentially expressed genes (DEGs), the RNA-Seq data were mapped to coding sequences with BWA MEM v. 0.7.12-r1039 [27], after summarizing mapped read counts with in-house Perl scripts, statistical tests were conducted with DESeq2 v. 1.8.2 [30]. Genes with BH method adjusted p-values 0.05 were defined as DEGs.

Query genes were obtained from UniProt KB [44], and were subjected to a BLASTP v2.2.22 [4] search against the amino acid sequences to determine tardigrade orthologs. For hormone receptor related gene identification, we screened the gene annotations for receptor genes that had Gene Ontology terms that contained “hormone”. Additionally, we obtained single embryo RNA-Seq data of H. exemplaris [26], and assembled a 3′ biased stranded transcriptome with Trinity v. 2.4.0 (default parameters). This transcriptome assembly was subjected to a BLASTn search against a previous assembly of the whole transcriptome of H. exemplaris, which then were further subjected to BLASTn searches against the coding sequences to identify the originating genes [48]. To quantify gene expression, the RNA-Seq data were mapped against the assembled 3′ biased stranded transcriptome with BWA MEM, and after summarizing mapped read counts with SAMtools idxstat v. 1.4 [28], TPM values were calculated in R.

Chemical exposure

Fenoxycarb (Wako, Japan) and 20-Hydroxyecdysone (20-E, Tokyo Chemical Industry Co., Ltd., Japan) were used as juvenile hormone analog and ecdysone, respectively. Chemicals were eluted in 15 μL 100% Ethanol (Nacalai tesque, Japan), and were further diluted with 45 ml of Volvic water to obtain the exposure concentrations (20-E: 8.32 μM, 100 μM, Fenoxycarb: 3.3 μM, 100 μM) following previous research on Daphnia [10, 42]. For controls, embryos were exposed to 15 μl 100% ethanol in 45 ml Volvic solution (0.03% Ethanol) at day 0. Embryos (0–3 days after oviposition) with exuvium (H. exemplaris) or freely laid eggs (R. varieornatus) were exposed to the chemical solution in 12 well dishes, and the number of days required for the eggs to hatch was observed for 6 days after oviposition. R. varieornatus embryos were exposed to only low concentrations. We also exposed embryos to 100 μM of ecdysone receptor antagonist Cucurbitacin B (Tokyo Chemical Industry Co., Ltd., Japan) to assess whether inhibiting the molting pathway affected embryo hatching [45, 50]. We assayed six and three clutches for low and high concentration exposures, respectively. Three clutches were assayed for the 100 μM Cucurbitacin B treatment. For R. varieornatus, we assayed three clutches for both chemicals. Treated samples (including control samples) were incubated at 18 °C. Furthermore, to validate if unhatched embryos were capable of hatching, embryos were washed with Volvic water three times to remove chemicals and then cultured until hatching could be confirmed. To assess whether the unhatched embryo had completed development, the cuticle shells of embryos that did not hatch at Egg 7d were peeled off with a 26G needle and were observed with light microscopy for claw formation.

Statistics

All statistics analysis were conducted in the R software. Standard deviation for each data are shown throughout this manuscript. Heatmaps and line plos were generated using R and ggplot2 library.

Availability of data and materials

All mRNA-Seq data have been uploaded to GEO under the accession ID GSE130111.

Abbreviations

20-E:

20-Hydroxyecdysone

ANOVA:

Analysis of Variance

EcR:

Ecdysone receptor

hpl:

hours post laying

JH:

Juvenile hormones

Kr-h1:

Krüppel homolog-1

SD:

Standard deviation

SEM:

Scanning Electron Microscope

SOM:

Self-organizing map

Tukey HSD:

Tukey Honestly Significant Difference

USP:

Ultraspiracle

References

  1. 1.

    Altiero T, Bertolani R, Rebecchi L. Hatching phenology and resting eggs in tardigrades. J Zool. 2010;280:290–6.

    Article  Google Scholar 

  2. 2.

    Altiero T, Rebecchi L, Bertolani R. Phenotypic variations in the life history of two clones of Macrobiotus richtersi (Eutardigrada, Macrobiotidae). Hydrobiologia. 2006;558:33–40.

    Article  Google Scholar 

  3. 3.

    Altiero T, Suzuki A, Rebecchi L. Reproduction, development and life cycles. In: Schill RO, editor. Ed. Water bears: the biology of Tardigrades. Basel: Springer International Publishing; 2019. p. 211–47.

    Google Scholar 

  4. 4.

    Altschul SF, Madden TL, Schaffer AA, Zhang J, Zhang Z, Miller W, Lipman DJ. Gapped BLAST and PSI-BLAST: a new generation of protein database search programs. Nucleic Acids Res. 1997;25:3389–402.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  5. 5.

    Andrews S. FastQC a quality-control tool for high-throughput sequence data; 2015.

    Google Scholar 

  6. 6.

    Arakawa K, Kido N, Oshita K, Tomita M. G-language genome analysis environment with REST and SOAP web service interfaces. Nucleic Acids Res. 2010;38:W700–5.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  7. 7.

    Arakawa K, Mori K, Ikeda K, Matsuzaki T, Kobayashi Y, Tomita M. G-language genome analysis environment: a workbench for nucleotide sequence data mining. Bioinformatics. 2003;19:305–6.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  8. 8.

    Arakawa K, Tomita M. G-language system as a platform for large-scale analysis of high-throughput omics data. J Pestic Sci. 2006;31:282–8.

    Article  CAS  Google Scholar 

  9. 9.

    Arakawa K, Yoshida Y, Tomita M. Genome sequencing of a single tardigrade Hypsibius dujardini individual. Sci Data. 2016;3:160063.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  10. 10.

    Baldwin WS, Bailey R, Long KE, Klaine S. Incomplete ecdysis is an indicator of ecdysteroid exposure in Daphnia magna. Environ Toxicol Chem. 2001;20:1564–9.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  11. 11.

    Bownes M, Rembold H. The titre of juvenile hormone during the pupal and adult stages of the life cycle of Drosophila melanogaster. Eur J Biochem. 1987;164:709–12.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  12. 12.

    Bray NL, Pimentel H, Melsted P, Pachter L. Near-optimal probabilistic RNA-seq quantification. Nat Biotechnol. 2016;34:525–7.

    Article  CAS  Google Scholar 

  13. 13.

    Gabriel WN, Goldstein B. Segmental expression of Pax3/7 and engrailed homologs in tardigrade development. Dev Genes Evol. 2007;217:421–33.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  14. 14.

    Gabriel WN, McNuff R, Patel SK, Gregory TR, Jeck WR, Jones CD, Goldstein B. The tardigrade Hypsibius dujardini, a new model for studying the evolution of development. Dev Biol. 2007;312:545–59.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  15. 15.

    Gordon R, Chippett J, Tilley J. Effects of two Carbamates on infective juveniles of Stemernema carpocapsae all strain and Steinernema feltiae Umeå strain. J Nematol. 1996;28:310–7.

    PubMed  PubMed Central  CAS  Google Scholar 

  16. 16.

    Gross V, Minich I, Mayer G. External morphogenesis of the tardigrade Hypsibius dujardini as revealed by scanning electron microscopy. J Morphol. 2017;278:563–73.

    PubMed  Article  PubMed Central  Google Scholar 

  17. 17.

    Hashimoto T, Horikawa DD, Saito Y, Kuwahara H, Kozuka-Hata H, Shin IT, Minakuchi Y, Ohishi K, Motoyama A, Aizu T, Enomoto A, Kondo K, Tanaka S, Hara Y, Koshikawa S, Sagara H, Miura T, Yokobori SI, Miyagawa K, Suzuki Y, Kubo T, Oyama M, Kohara Y, Fujiyama A, Arakawa K, Katayama T, Toyoda A, Kunieda T. Extremotolerant tardigrade genome and improved radiotolerance of human cultured cells by tardigrade-unique protein. Nat Commun. 2016;7:12808.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  18. 18.

    Horikawa DD, Kunieda T, Abe W, Watanabe M, Nakahara Y, Yukuhiro F, Sakashita T, Hamada N, Wada S, Funayama T, Katagiri C, Kobayashi Y, Higashi S, Okuda T. Establishment of a rearing system of the extremotolerant tardigrade Ramazzottius varieornatus: a new model animal for astrobiology. Astrobiology. 2008;8:549–56.

    PubMed  Article  CAS  Google Scholar 

  19. 19.

    Ito M, Saigo T, Abe W, Kubo T, Kunieda T. Establishment of an isogenic strain of the desiccation-sensitive tardigrade Isohypsibius myrops (Parachela, Eutardigrada) and its life history traits. Zool J Linnean Soc. 2016;178:863–70.

    Article  Google Scholar 

  20. 20.

    Karlson P. On the hormonal control of insect metamorphosis: a historical review. Int J Dev Biol. 1996;40:93–6.

    PubMed  CAS  Google Scholar 

  21. 21.

    Kato Y, Kobayashi K, Oda S, Tatarazako N, Watanabe H, Iguchi T. Cloning and characterization of the ecdysone receptor and ultraspiracle protein from the water flea Daphnia magna. J Endocrinol. 2007;193:183–94.

    PubMed  Article  CAS  Google Scholar 

  22. 22.

    Kayukawa T, Jouraku A, Ito Y, Shinoda T. Molecular mechanism underlying juvenile hormone-mediated repression of precocious larval-adult metamorphosis. Proc Natl Acad Sci U S A. 2017;114(5):1057–62.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  23. 23.

    Keilin D. The Leeuwenhoek lecture - the problem of Anabiosis or latent life - history and current concept. Proc Royal Soc B Biol Sci. 1959;150:149–91.

    Article  CAS  Google Scholar 

  24. 24.

    Koziol U. Precursors of neuropeptides and peptide hormones in the genomes of tardigrades. Gen Comp Endocrinol. 2018;267:116–27.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  25. 25.

    Langmead B, Salzberg SL. Fast gapped-read alignment with bowtie 2. Nat Methods. 2012;9:357–9.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  26. 26.

    Levin M, Anavy L, Cole AG, Winter E, Mostov N, Khair S, Senderovich N, Kovalev E, Silver DH, Feder M, Fernandez-Valverde SL, Nakanishi N, Simmons D, Simakov O, Larsson T, Liu SY, Jerafi-Vider A, Yaniv K, Ryan JF, Martindale MQ, Rink JC, Arendt D, Degnan SM, Degnan BM, Hashimshony T, Yanai I. The mid-developmental transition and the evolution of animal body plans. Nature. 2016;531:637–41.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  27. 27.

    Li H, Durbin R. Fast and accurate short read alignment with burrows-wheeler transform. Bioinformatics. 2009;25:1754–60.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  28. 28.

    Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R, Genome Project Data Processing, S. The sequence alignment/map format and SAMtools. Bioinformatics. 2009;25:2078–9.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  29. 29.

    Li K, Jia QQ, Li S. Juvenile hormone signaling - a mini review. Insect Sci. 2019;26(4):600–6.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  30. 30.

    Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  31. 31.

    Metwally HMS, Abdel-Hakim E, Abdou WL. Effects of the Juvenoid (Fenoxycarb) on Entomopathogenic nematode Steinernema carpocapsae S2 strain. Biosci Res. 2017;14:966–75.

    Google Scholar 

  32. 32.

    Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26:841–2.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  33. 33.

    Ray R, Pandey P. piRNA analysis framework from small RNA-Seq data by a novel cluster prediction tool - PILFER. Genomics. 2018;110:355–65.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  34. 34.

    Rosenkranz D, Zischler H. proTRAC--a software for probabilistic piRNA cluster detection, visualization and analysis. BMC Bioinformatics. 2012;13:5.

    PubMed  PubMed Central  Article  Google Scholar 

  35. 35.

    Sarkies P, Selkirk ME, Jones JT, Blok V, Boothby T, Goldstein B, Hanelt B, Ardila-Garcia A, Fast NM, Schiffer PM, Kraus C, Taylor MJ, Koutsovoulos G, Blaxter ML, Miska EA. Ancient and novel small RNA pathways compensate for the loss of piRNAs in multiple independent nematode lineages. PLoS Biol. 2015;13:e1002061.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  36. 36.

    Schill RO, Fritz GB. Desiccation tolerance in embryonic stages of the tardigrade. J Zool. 2008;276:103–7.

    Article  Google Scholar 

  37. 37.

    Schumann I, Kenny N, Hui J, Hering L, Mayer G. Halloween genes in panarthropods and the evolution of the early moulting pathway in Ecdysozoa. R Soc Open Sci. 2018;5:180888.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  38. 38.

    Sekimoto T, Iwami M, Sakurai S. Coordinate responses of transcription factors to ecdysone during programmed cell death in the anterior silk gland of the silkworm, Bombyx mori. Insect Mol Biol. 2006;15:281–92.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  39. 39.

    Smith FW, Boothby TC, Giovannini I, Rebecchi L, Jockusch EL, Goldstein B. The compact body plan of Tardigrades evolved by the loss of a large body region. Curr Biol. 2016;26:224–9.

    PubMed  Article  CAS  Google Scholar 

  40. 40.

    Sonobe H, Yamada R. Ecdysteroids during early embryonic development in silkworm Bombyx mori: metabolism and functions. Zool Sci. 2004;21:503–16.

    PubMed  Article  CAS  Google Scholar 

  41. 41.

    Tanaka S, Tanaka J, Miwa Y, Horikawa DD, Katayama T, Arakawa K, Toyoda A, Kubo T, Kunieda T. Novel mitochondria-targeted heat-soluble proteins identified in the anhydrobiotic Tardigrade improve osmotic tolerance of human cells. PLoS One. 2015;10:e0118272.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  42. 42.

    Tatarazako N, Oda S, Watanabe H, Morita M, Iguchi T. Juvenile hormone agonists affect the occurrence of male Daphnia. Chemosphere. 2003;53:827–33.

    PubMed  Article  CAS  Google Scholar 

  43. 43.

    Tenlen JR, McCaskill S, Goldstein B. RNA interference can be used to disrupt gene function in tardigrades. Dev Genes Evol. 2013;223:171–81.

    PubMed  Article  CAS  Google Scholar 

  44. 44.

    UniProt C. UniProt: a worldwide hub of protein knowledge. Nucleic Acids Res. 2019;47:D506–15.

    Article  CAS  Google Scholar 

  45. 45.

    Vellichirammal NN, Gupta P, Hall TA, Brisson JA. Ecdysone signaling underlies the pea aphid transgenerational wing polyphenism. Proc Natl Acad Sci U S A. 2017;114:1419–23.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  46. 46.

    Wang J, Zhang P, Lu Y, Li Y, Zheng Y, Kan Y, Chen R, He S. piRBase: a comprehensive database of piRNA sequences. Nucleic Acids Res. 2019;47:D175–80.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

  47. 47.

    Yamaguchi A, Tanaka S, Yamaguchi S, Kuwahara H, Takamura C, Imajoh-Ohmi S, Horikawa DD, Toyoda A, Katayama T, Arakawa K, Fujiyama A, Kubo T, Kunieda T. Two novel heat-soluble protein families abundantly expressed in an anhydrobiotic tardigrade. PLoS One. 2012;7:e44209.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  48. 48.

    Yoshida Y, Koutsovoulos G, Laetsch DR, Stevens L, Kumar S, Horikawa DD, Ishino K, Komine S, Kunieda T, Tomita M, Blaxter M, Arakawa K. Comparative genomics of the tardigrades Hypsibius dujardini and Ramazzottius varieornatus. PLoS Biol. 2017;15:e2002266.

    PubMed  PubMed Central  Article  CAS  Google Scholar 

  49. 49.

    Yoshida Y, Sugiura K, Tomita M, Matsumoto M, Arakawa K. Data from “Comparison of the transcriptomes of two tardigrades with different hatching coordination”: Figshare; 2019. https://doi.org/10.6084/m9.figshare.10163933.v1. Deposited 2 Nov.

  50. 50.

    Zou C, Liu G, Liu S, Liu S, Song Q, Wang J, Feng Q, Su Y, Li S. Cucurbitacin B acts a potential insect growth regulator by antagonizing 20-hydroxyecdysone activity. Pest Manag Sci. 2018;74:1394–403.

    PubMed  Article  CAS  PubMed Central  Google Scholar 

Download references

Acknowledgements

We thank Nozomi Abe, Yuki Takai, and Naoko Ishii for tardigrade sample preparation and transcriptome sequencing. We also thank Hiroki Minato and Ryusei Otsuka for suggesting the optimal concentration of hormones exposure. We also thank the anonymous reviewers for their critical comments. We are grateful to Dr. James Fleming for help in improving the English in this manuscript.

Funding

YY and KS are funded by Japan Society for the Promotion of Science (JSPS), JP18J21155 for YY, JP18J21345 for KS. This work is supported by KAKENHI Grant-in-Aid for Scientific Research (B) from the JSPS (grant no. 17H03620, for KA), and partly by research funds from the Yamagata Prefectural Government and Tsuruoka City, Japan for KA and MT. The funders had no role in the design of the study and collection, analysis, and interpretation of data and in writing the manuscript.

Author information

Affiliations

Authors

Contributions

KA conceived the project and performed RNA sequencing, and YY performed all bioinformatics analysis. MM and KS designed the life history and chemical exposure experiments, and KS conducted these experiments. MT managed the compute resources and provided supervision for YY. YY and KS drafted the manuscript, and all authors have read and approved the manuscript.

Corresponding author

Correspondence to Kazuharu Arakawa.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Additional file 1:

Figure S1. Body length in developing individuals of H. exemplaris. The body length was quantified for new hatchlings and observed for 28 days.

Additional file 2:

Figure S2. HOX genes were upregulated around Egg 2-3d in both species. Z-scaled TPM values of HOX genes identified in Yoshida et al. [48] were visualized as a heatmap. E: Egg, B: Juvenile, Active/Tun: Adult stages.

Additional file 3:

Figure S3. Expression of HOX genes in the CEL-Seq data set. CEL-Seq reads were assembled to construct a 3′ strand biased transcriptome, used to identify HOX genes. The CEL-Seq reads were mapped against this assembly to calculate gene expression (TPM).

Additional file 4:

Figure S4. Expression profiles of CAHS and SAHS orthologs. TPM values of (A) CAHS and (B) SAHS genes were plotted with standard deviation as error bars.

Additional file 5:

Figure S5. Profiles of anti-oxidative stress related genes. Z-scaled TPM values of catalase (CATE), glutathione S-transferase (GST), superoxide dismutase (SOD), and Thioredoxin reductase (TRX) of both species were plotted as a heatmap.

Additional file 6:

Figure S6. Duration of H. exemplaris embryo exposed to low concentration chemicals. Days required for hatching in embryo exposed to low concentrations. NH: Not hatched. Error bars indicates standard variation.

Additional file 7:

Figure S7. Exposure to low concentration chemicals in both tardigrades. Embryo were exposed to low concentrations of Fenoxycarb and 20-E at the indicated days, and the hatching ratio was recorded. Slightly higher variation was observed in the H. exemplaris embryos, but not in R. varieronatus.

Additional file 8:

Figure S8. Duration of H. exemplaris embryogenesis exposed to high-concentration chemicals. Days required for hatching in embryo exposed to high concentrations. NH: Not hatched. Error bars indicates standard variation.

Additional file 9:

Table S1. Summary of RNA-Seq data.

Additional file 10.

Expression profiles of genes in group 1 and 2 of SOM clustering in H. exemplaris. https://github.com/abs-yy/Yoshida_and_Sugiura_etal_2018/blob/master/Sdata_1_hypsibius_som_group1ANND2.xlsx.

Additional file 11.

Excel file containing the expression profiles of genes focused in this analysis. https://github.com/abs-yy/Yoshida_and_Sugiura_etal_2018/blob/master/Sdata_2_geneexpression.xlsx.

Additional file 12.

Annotations of genes overlapping piRNA loci predicted in both PILFER and proTRAC. https://github.com/abs-yy/Yoshida_and_Sugiura_etal_2018/blob/master/piRNA_clusters_intersect_geneAnnotations.txt.

Additional file 14.

Additional text on the analysis pf piRNA genes and anhydrobiosis related genes.

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Yoshida, Y., Sugiura, K., Tomita, M. et al. Comparison of the transcriptomes of two tardigrades with different hatching coordination. BMC Dev Biol 19, 24 (2019). https://doi.org/10.1186/s12861-019-0205-9

Download citation

Keywords

  • Tardigrades
  • Developmental stages
  • RNA-Seq
  • Ecdysone