Skip to main content

The developmental transcriptomes of two sea biscuit species with differing larval types



Larval developmental patterns are extremely varied both between and within phyla, however the genetic mechanisms leading to this diversification are poorly understood. We assembled and compared the developmental transcriptomes for two sea biscuit species (Echinodermata: Echinoidea) with differing patterns of larval development, to provide a resource for investigating the evolution of alternate life cycles. One species (Clypeaster subdepressus) develops via an obligately feeding larva which metamorphoses 3–4 weeks after fertilization; the other (Clypeaster rosaceus) develops via a rare, intermediate larval type—facultative feeding— and can develop through metamorphosis entirely based on egg provisioning in under one week.


Overall, the two transcriptomes are highly similar, containing largely orthologous contigs with similar functional annotation. However, we found distinct differences in gene expression patterns between the two species. Larvae from C. rosaceus, the facultative planktotroph, turned genes on at earlier stages and had less differentiation in gene expression between larval stages, whereas, C. subdepressus showed a higher degree of stage-specific gene expression.


This study is the first genetic analysis of a species with facultatively feeding larvae. Our results are consistent with known developmental differences between the larval types and raise the question of whether earlier onset of developmental genes is a key step in the evolution of a reduced larval period. By publishing a transcriptome for this rare, intermediate, larval type, this study adds developmental breadth to the current genetic resources, which will provide a valuable tool for future research on echinoderm development as well as studies on the evolution of development in general.


Over the last century, echinoid echinoderms have come to represent an iconic model system for studying questions about the evolution of development, due in large part to their diverse larval forms, even among closely related species [1,2,3,4,5,6,7,8,9,10,11,12,13]. However, only recently has it become possible to decipher the genetic changes underlying evolutionary shifts in developmental mode [14,15,16,17,18]. Here, we present the developmental transcriptomes of two closely related echinoids as a foundation for analyzing how patterns of gene expression change during the course of larval evolution.

Despite their great morphological diversity, marine larvae are typically classified into two groups: planktotrophic larvae, that develop from small eggs and must feed in order to develop through metamorphosis, and lecithotrophic larvae, that do not feed as larvae but instead develop entirely based on energy provided in their large yolky eggs [7]. Feeding, planktotrophic larvae presumably represent the ancestral state for all echinoderms, whereas non-feeding larvae have evolved convergently dozens of times [13, 19]. While the evolution of non-feeding larvae from feeding larvae involves many morphological and functional changes, including (1) an increase in egg size; (2) reduction in developmental time; (3) accelerated development of a juvenile body; and (4) the reduction—or even complete loss—of larval feeding structures, the transition can occur rapidly [20,21,22,23,24,25,26,27].

There are a handful of intermediate larval types between these two extreme modes, which have proven invaluable for understanding the speed, manner, and order in which developmental changes occur during the evolution of non-feeding forms [26, 28,29,30]. The sea biscuit Clypeaster rosaceus is perhaps the best-studied echinoderm with an intermediate larval type [27, 31,32,33]. Like obligate planktotrophs, including its closely related congener C. subdepressus, C. rosaceus larvae develop through a pluteus larval form with functioning feeding structures [31]. However, unlike obligately planktotrophic larvae, C. rosaceus larvae develop from large eggs, can complete metamorphosis without feeding [33] and develop in under a week, as opposed to the three weeks C. subdepressus takes [31, 34]. Furthermore, C. rosaceus larvae begin forming structures that develop into a juvenile shortly after gastrulation, like species with non-feeding development, but unlike species with obligately feeding larvae that typically do not produce juvenile structures until late larval stages [26]. Lastly, the larval feeding structures of C. rosaceus are smaller and less efficient relative to obligate planktotrophs [31, 35]. Consequently, studies on the development of C. rosaceus, as well as other species with intermediate larval types, fill a void in our understanding of the order and timing of morphological and genetic changes during the evolutionary transition from feeding to non-feeding larvae [28, 29].

Published transcriptomes for echinoderm species with both planktotrophic [16, 17, 36, 37] and lecithotrophic [38] development reveal dramatic changes in gene expression and gene regulatory networks between the two [17]. However, there are currently no assembled transcriptomes for species with intermediate larval types. We therefore assembled transcriptomes for C. rosaceus, along with a closely related congener, C. subdepressus [33]. These two species are an ideal study system for understanding the order and manner in which divergent patterns of gene expression accumulate during the evolutionary transition between larval types not only because of their close phylogenetic relationship, but also due to the extensive previous work on their early development [31, 32, 35], and the fact that they can produce viable hybrid larvae [34].


Tissue collection

We collected adult sea biscuits (C. rosaceus and C. subdepressus) from subtidal sand beds at Bocas Del Toro, Panama during the months of August and September 2013, and maintained them in ambient conditions in outdoor flow-through tanks at the Bocas Del Toro Research Station of the Smithsonian Tropical Research Institute. To generate embryos and larvae, we established three conspecific replicate crosses for each species. In each replicate cross, we fertilized eggs from a single female with sperm from a single male in 100 ml finger bowls, so that the larvae of each species were full siblings within a replicate, but were unrelated across replicates.

We reared larvae according to the methods described in Armstrong and Lessios [34]. We set up 1000 ml tri-pour beakers with sea-water filtered through a 0.45μm membrane. After embryos hatched, we added 200 larvae of either C. rosaceus or C. subdepressus to each beaker. We maintained larvae on a Strathmann stirring rack [8] at ambient temperature conditions. Every other day, we changed 80% of the water in each beaker and fed larvae a mix of the two algal species Isochrysis galbana and Rhodomonas spp. at a concentration of 5000 cells/ml.

We collected RNA samples at four developmental stages: unfertilized eggs (0 h post-fertilization: 0 hpf), gastrulae (20 hpf for both species), four-armed larvae (2 days post-fertilization (dpf) in C. rosaceus; 3 dpf in C. subdepressus), and 8-armed larvae (4 dpf in C. rosaceus; 10 dpf in C. subdepressus) [developmental staging from Armstrong and Lessios [34]]. For each RNA sample, we preserved 100 individuals (eggs, gastrulae, or larvae) in RNAlater. We stored the samples at −20 °C while at Bocas Del Toro, followed by −80 °C for long-term storage at UC Davis.

RNA extraction and sequencing

We extracted total RNA using Trizol, followed by a DNAse treatment using Ambion’s Turbo DNA-free kit. We quantified RNA using a Nanodrop and assessed RNA quality on an Agilent 2100 Bioanalyzer; we discarded any RNA samples with a RNA integrity number (RIN) < 7. From the resulting RNA samples, we created six replicate RNAseq libraries (two libraries for each of the three replicate experiments) for each developmental stage in each species, following standard protocols with the NEB Ultra Directional kit. We multiplexed the libraries together by using NEB multiplex primers and sequenced the libraries across three lanes of Paired-End100 on an Illumina HiSeq2500 machine at the UC Davis Genome Center.

Assembly and annotation

We removed adapter sequences and low-quality reads from the raw Illumina reads with the program Trimmomatic (v0.36) using default parameters; we retained sequences that were at least 36 bp after trimming. We combined all of the remaining reads for each species and assembled transcriptomes using Trinity (v2.4.0) set to default parameters with a minimum contig length of 300 bp [39]. To account for redundant transcripts due to within-species polymorphisms or sequencing errors, we clustered the resultant Trinity outputs at 97% similarity using the program CD-hit (v4.6.5) [40]. For each clustered transcriptome, we identified contigs with coding potential through the program TransDecoder [41]. Finally, we used the script from the Transcriptome Utilities package to remove contigs originating from bacterial and viral contaminants [42].

We used the programs CEGMA (Core Eukaryotic genes mapping approach) [43] and BUSCO (Benchmarking universal single-copy orthologs) [44] to assess the completeness of core eukaryotic and metazoan proteins, respectively, in each transcriptome. To identify orthologous contigs between the transcriptomes, we used the program OrthoFinder [45] with default parameters and compared our two assemblies along with the purple urchin transcriptome, Strongylocentrotus purpuratus (downloaded from ) [46]. We functionally annotated the transcriptomes using the program Blast2GO (v4.1) [47] with the UniProt BLAST database. We ran the program Interpro scan to obtain functional information, including protein families, Panther biological process categories [48], and gene ontology (GO) terms [49]. We compared the overall GO distribution between the two transcriptomes using the R package GOstats to test for any significantly enriched biological process GO terms [50].

Gene expression analysis

For differential expression analysis, we only retained libraries with > 3 million reads after quality filtering. We mapped our filtered reads back to their respective transcriptome using the program RSEM [41] to obtain expression tables. We normalized the data with a TMM (trimmed mean of m-values) method [51]. Using the normalized count data, we analyzed the Euclidean distances between the libraries from each species using a heatmap in DESeq [52]; we also set up a generalized linear model in EdgeR [53] with developmental stage and replicate experiment as fixed effects. From this model, we identified genes that were significantly differentially expressed in each species between developmental stages, with a false discovery rate of 0.05. Using the resulting lists of up-regulated genes at each stage, we identified any biological process GO terms that were over expressed between developmental time points. We used a hypergeometric test in the R package GOstats for this analysis [50].



Using three lanes of HiSeq PE100, we generated 240 and 260 million reads for C. rosaceus and C. subdepressus, respectively, with a mean phred quality score > 36. After removing low-quality reads and adapters in trimmomatic, we retained > 95% (237 and 252 million reads, respectively) to assemble each transcriptome. The final assemblies had 33,190 and 36,059 contigs for C. rosaceus and C. subdepressus, with N50’s of 2769 bp and 2430 bp, respectively (Table 1). Both transcriptomes contained complete orthologs for between 96 and 99% of core eukaryotic and metazoan proteins when analyzed in the programs CEGMA and BUSCO respectively (Table 1). We identified orthologs between the two transcriptomes, along with the purple sea urchin transcriptome, using OrthoFinder (Emms and Kelly 2015). Orthofinder generated 15,038 orthogroups with 10,779 (72%) containing all three species (Fig. 1). The two Clypeaster species shared an additional 2625 orthogroups exclusively with each other, more than triple the number of orthogroups either species exclusively shared with S. purpuratus (Fig. 1).

Table 1 Assembly statistics for final transcriptomes of C. rosaceus and C. subdepressus
Fig. 1
figure 1

Summary of orthogroups between Clypeaster rosaceus (red), C. subdepressus (yellow), and Strongylocentrotus purpuratus (purple). Orthogroups with all three species present are shown overlapping in the middle (10,779). Orthogroups shared exclusively by the Clypeaster species are shown in orange (2625)

Our BLAST search yielded at least one significant hit for 25,440 (77%) C. rosaceus contigs and 28,947 (80%) C. subdepressus contigs. Additionally, we obtained at least one GO term for 16,240 (49%) C. rosaceus contigs and 17,278 (48%) C. subdepressus contigs. Our analysis revealed no significant enrichment for any biological process GO terms between the two transcriptomes (Fig. 2). In both transcriptomes, cellular process, metabolic process, and response to stimulus were the largest GO terms, respectively (Fig. 2).

Fig. 2
figure 2

Distribution of biological process gene ontology terms in each transcriptome. The category and percent of GO terms belonging to that category are labeled for both Clypeaster rosaceus (a) and C. subdepressus (b)

Differential gene expression

After quality filtering and discarding any libraries with < 3 million reads, we retained 20 libraries per species, consisting of at least four replicates per developmental stage in each species. We used these libraries to analyze differential expression across developmental stages within each species.

To determine which libraries were most similar, we analyzed the Euclidean distance between the samples from each species, using the heatmap function in DESeq [52]. From this analysis it was evident that, in both species, the eggs were the most distinct developmental time point from all other samples, followed by the gastrulae, and subsequently the two larval stages clustered together (Fig. 3). In C. subdepressus, the larval stages clearly clustered into two separate groups based on developmental stage. However, in the C. rosaceus samples, both larval stages were intermixed without a clear pattern. In both species, there was an effect of experiment, as replicates from the same experiment generally clustered together.

Fig. 3
figure 3

Heatmaps of the Euclidean distance between samples for each Clypeaster rosaceus (a), and C. subdepressus (b). Dark purple colors represent highly similar samples, while lighter, yellow colors reflect more distant samples. The order of samples is denoted below each column with the sample names indicating first the developmental stage, followed by experiment (1–3), and replicate within experiment (a or b). The branching pattern at the top of each plot represents relative distances between samples

Next, we identified differentially expressed genes (DEGs) using a generalized linear model in EdgeR, with both developmental stage and replicate experiment as fixed effects (Fig. 4). In both species, the greatest differences in expression occurred between the egg and gastrula stages, followed by the gastrula to 4-armed larva; and the smallest amount of differential expression occurred between the 4-arm and 8-arm larval stages (Fig. 4). More genes were differentially expressed during the egg-to-gastrula transition in C. rosaceus relative to C. subdepressus. In C. rosaceus 39.4% (13,062) of genes were expressed at statically different levels between the egg and gastrula stage as opposed to 34.6% (12,490) in C. subdepressus. However, the opposite was true in the other two contrasts (gastrula to 4-armed larva, and 4-armed to 8-armed larva) where there were more DEGs in C. subdepressus relative to C. rosaceus. This was most striking in the comparison between the four and eight-armed larvae, where only seven (0.02%) genes were differentially expressed in C. rosaceus whereas 633 genes (1.8%) were in C. subdepressus. There was a significant effect of replicate experiment, with 4.7% of C. rosaceus genes and 4.8% of C. subdepressus genes showing differences in gene expression between our three replicates.

Fig. 4
figure 4

Scatterplots of stage specific gene expression. a-c Clypeaster rosaceus, (d-f) C. subdepressus. Expression levels are displayed in counts per million contrasted pairwise between developmental stages as labeled on the x- and y-axes. Red points represent genes identified as significantly differentially expressed, using a false discovery rate of 0.05, between stages while black dots are not significantly different. The number (and percent) of differentially expressed genes is indicated in each plot

We used the R-package GOstats to identify GO terms that were overrepresented in each list of differentially expresses genes at each developmental stage. This analysis revealed more significantly overrepresented GO terms at every developmental stage of C. subdepressus than in C. rosaceus (Additional file 1: Table S1). There was a high amount of overlap in the enriched GO categories between species at both the egg and gastrula stages; over 90% of significantly enriched GO terms in the eggs or gastrulae of C. rosaceus were also significantly enriched at the same stages of C. subdepressus. This trend, however, did not hold true when comparing significant GO terms in the larval stages (Additional file 1: Table S1).

In the eggs of both species, the GO terms with the greatest overrepresentation were related to protein modification, cell communication, signaling, and regulation of biological processes (GO:0036211, GO:0050794, GO:0051716, GO:0050789, GO:0007154, GO:0023052). In the gastrulae of both species the GO terms most enriched were related to translation, gene expression, biosynthetic processes (including ribosome biogenesis), mRNA processing, and protein folding (GO:0006412, GO:0010467, GO:0009058, GO:0022613, GO:0042254, GO:0006397). At the 4-arm larval stage, the only GO term enriched in both species was for transmembrane transport (GO:0055085). Additionally, cell communication, signaling, and cell cycle were enriched in C. subdepressus 4-armed larvae (GO:0007154, GO:0023052; GO:0007049, GO:0000278), whereas in C. rosaceus biosynthetic and metabolic processes were enriched (GO:0005975, GO:0044249, GO:0043043). At the 8-arm larval stage of C. subdepressus, metabolic processes were statistically enriched (GO:0044237, GO:0008152, GO:0006629), whereas no GO terms were significantly enriched in the 8-arm stage of C. rosaceus (for a full list of significant GO terms see Additional file 1: Table S1).


The developmental transcriptomes we assembled for C. rosaceus and C. subdepressus are highly consistent with one another in terms of assembly statistics, functional annotation and completeness of core eukaryotic proteins. The majority of contigs had orthologs in the other transcriptome, including many shared orthogroups that were absent from the Strongylocentrotus purpuratus outgroup. The two Clypeaster species diverged about 8 million years ago [33], which is consistent with the low level of differentiation between their transcriptomes.

Both transcriptomes showed the largest differences between ontological states in gene expression at the earliest developmental stages, with diminishing numbers of DEGs as development progressed. This result is expected, as a vast amount of developmental and morphological rewiring occurs between the egg and gastrula stage, compared to the minimal differences between larval stages [38, 54]. In echinoids, zygotic transcription begins as early as the first cleavage [55,56,57], so that while the transcripts in unfertilized eggs are entirely maternally deposited, the transcripts in the gastrulae should be entirely zygotic. Additionally, many studies report a spike of transcripts involved in ribosome production, translation, and protein synthesis in early embryonic phases, as we observed in the gastrulae of both species [36, 57]. Overall, these patterns were consistent between our species of Clypeaster, and with studies of embryonic development in other species [18, 36, 37].

Obligately planktotrophic larvae, which must acquire food to develop, focus their energy on feeding and producing feeding structures during the four-arm larval stage, while postponing the production of a juvenile rudiment [26, 58]. Alternatively, the facultatively feeding larvae of C. rosaceus begin developing juvenile structures shortly after gastrulation rather than solely devoting their energy towards feeding at early larval stages [26]. Due to these developmental differences, gene expression should be more dynamic in C. subdepressus where the embryonic and larval stages have more distinct functions. Not surprisingly, the most striking difference we found between these two species was that Clypeaster subdepressus displayed more dynamic, stage-specific, gene expression patterns—as indicated by our differential expression results, cluster analysis, and GO analysis—than C. rosaceus, particularly during the later larval stages. The only time point where we observed a higher turnover of differentially expressed genes in C. rosaceus relative to C. subdepressus was between the egg and gastrula stages; consistent with the fact that facultatively feeding larvae initiate the development of a juvenile body at the gastrula stage, whereas obligately planktotrophic larvae do not begin this process until late larval stages. Thus, just as accelerated developmental processes are one of the first steps in the evolution of non-feeding larvae [26, 27], our data suggest that accelerated gene activation is a key initial step in the loss of larval feeding.


The transcriptomes we present add both developmental and phylogenetic breadth to the currently available larval transcriptomes. Both assemblies were highly complete and similar to one another but displayed distinct expression patterns that fit with our understanding of the developmental differences between these two species. Our results indicate that similar genes are expressed throughout development in C. rosaceus and C. subdepressus. However, the timing of gene expression seems accelerated in the intermediate, facultatively feeding larvae, of C. rosaceus relative to the obligately feeding, C. subdepressus larvae. Future studies using additional species with intermediate larval types, such as the heart urchin Brisaster latifrons [11], could further explore whether these differences are indeed related to life history rather than species-specific differences. However, as C. rosaceus and B. latifrons are the only with only two documented echinoids with facultatively feeding larval types, such a comparative analysis would be limited in scope. Alternatively, it is possible to hybridize C. rosaceus and C. subdepressus. Each reciprocal hybrid cross has a similar zygotic genome, with half of their DNA coming from each species, regardless of which species the egg or sperm came from. However, offspring developmental mode depends on which species provides the egg [34]. Thus, the hybrid offspring provide further contrasts of developmental type while accounting for species-specific genomic differences. As such, future studies using these hybrid offspring provide a unique opportunity to reveal critical insights into the patterns and regulatory mechanisms involved in larval evolution.



Basic local alignment search tool


base pair


Benchmarking universal single-copy orthologs


Core Eukaryotic genes mapping approach


Differentially expressed gene


Deoxynucleic acid


days post-fertilization


Gene ontology


hours post-fertilization


RNA integrity number


Ribonucleic acid


Trimmed mean of m-values


  1. Vernon HM. The effect of environment on the development of echinoderm larvae: an experimental inquiry into the causes of variation. Philos Trans R Soc Lond B. 1895;186:577–632.

    Article  Google Scholar 

  2. Mortensen T. On the development of some British echinoderms. J Mar Biol Assoc U K. 1913;10:1–18.

    Article  Google Scholar 

  3. Just EE. The fertilization reaction in Echinarachnius parma: I. Cortical response of the egg to insemination. Biol Bull. 1919;36:1–10.

    Article  CAS  Google Scholar 

  4. Thorson G. Reproductive and larval ecology of marine bottom invertebrates. Biol Rev. 1950;25:1–45.

    Article  PubMed  CAS  Google Scholar 

  5. Okazaki K, Dan K. The metamorphosis of partial larvae of Peronella japonica Mortensen, a sand dollar. Biol Bull. 1954;106:83–99.

    Article  Google Scholar 

  6. Strathmann RR. The feeding behavior of planktotrophic echinoderm larvae: mechanisms, regulation, and rates of suspension feeding. J Exp Mar Biol Ecol. 1971;6:109–60.

    Article  Google Scholar 

  7. Strathmann RR. The evolution and loss of feeding larval stages of marine invertebrates. Evolution. 1978;32:894–906.

    Article  PubMed  Google Scholar 

  8. Strathmann MF. Reproduction and development of marine invertebrates of the northern pacific coast: data and methods for the study of eggs, embryos, and larvae: University of Washington Press; 1987.

  9. Christiansen FB, Fenchel TM. Evolution of marine invertebrate reproductive patterns. Theor Popul Biol. 1979;16:267–82.

    Article  PubMed  CAS  Google Scholar 

  10. Wray GA, Raff RA. The evolution of developmental strategy in marine invertebrates. Trends in Ecology and Evolution. 1991;6:45–50.

    Article  PubMed  CAS  Google Scholar 

  11. Hart MW. Evolutionary loss of larval feeding: development, form and function in a facultatively feeding larva, Brissaster latifrons. Evolution. 1996;50:174–87.

    Article  PubMed  Google Scholar 

  12. McEdward LR. Reproductive strategies of marine benthic invertebrates revisited: facultative feeding by planktotrophic larvae. Am Nat. 1997;150:48–72.

    Article  PubMed  CAS  Google Scholar 

  13. Jeffery CH, Emlet RB. Macroevolutionary consequences of developmental mode in temnopleurid echinoids from the tertiary of southern Australia. Evolution. 2003;57:1031–48.

    Article  PubMed  Google Scholar 

  14. Nielsen MG, Wilson KA, Raff EC, Raff RA. Novel gene expression patterns in hybrid embryos between species with different modes of development. Evol Dev. 2000;2:133–44.

    Article  PubMed  CAS  Google Scholar 

  15. McIntyre DC, Lyons DC, Martik M, McClay DR. Branching out: origins of the sea urchin larval skeleton in development and evolution. Genesis. 2014;52:173–85.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Tu Q, Cameron RA, Davidson EH. Quantitative developmental transcriptomes of the sea urchin. Dev Biol. 2014;385:160–7.

    Article  PubMed  CAS  Google Scholar 

  17. Israel JW, Martik ML, Byrne M, Raff EC, Raff RA, McClay DR, Wray GA. Comparative developmental transcriptomics reveals rewiring of a highly conserved gene regulatory network during a major life history switch in the sea urchin genus Heliocidaris. PLoS Biol. 2016;14:e1002391.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  18. Cary GA, Hinman VF. Echinoderm development and evolution in the post-genomic era. Dev Biol. 2017;427:203–11.

    Article  PubMed  CAS  Google Scholar 

  19. Raff EC, Popodi EM, Kauffman JS, Sly BJ, Turner FR, Morris VB, Raff RA. Regulatory punctuated equilibrium and convergence in the evolution of developmental pathways in direct developing sea urchins. Evol Dev. 2003;5:478–93.

    Article  PubMed  Google Scholar 

  20. Wray GA. Parallel evolution of nonfeeding larvae in echinoids. Syst Biol. 1996;45:308–22.

    Article  Google Scholar 

  21. Zigler KS, Raff EC, Popodi E, Raff RA, Lessios HA. Adaptive evolution of bindin in the genus Heliocidaris is correlated with the shift to direct development. Evolution. 2003;57:2293–302.

    Article  PubMed  CAS  Google Scholar 

  22. Raff RA, Byrne M. The active evolutionary lives of echinoderm larvae. Heredity. 2006;97:244–52.

    Article  PubMed  CAS  Google Scholar 

  23. Puritz JB, Keever CC, Addison JA, Byrne M, Hart MW, Grosberg RK, Toonen RJ. Extraordinarily rapid life-history divergence between Cryptasterina Sea star species. Proc R Soc Lond B Biol Sci. 2012;279:3914–22.

    Article  Google Scholar 

  24. Henry JJ, Klueg KM, Raff RA. Evolutionary dissociation between cleavage, cell lineage and embryonic axes in sea urchin embryos. Development. 1992;114:931–8.

    PubMed  CAS  Google Scholar 

  25. Raff EC, Popodi EM, Sly BJ, Turner R, Vilinski JT, Raff RA. A novel ontogenetic pathway in hybrid embryos between species with different modes of development. Development. 1999;126:1937–45.

    PubMed  CAS  Google Scholar 

  26. Snoke-Smith MS, Zigler KS, Raff RA. Evolution of direct-developing larvae: selection vs. loss. BioEssays. 2007;29:566–71.

    Article  Google Scholar 

  27. Zigler KS, Raff RA. A shift in germ layer allocation is correlated with large egg size and facultative planktotrophy in the echinoid Clypeaster rosaceus. Biol Bull. 2013;224:192–9.

    Article  PubMed  Google Scholar 

  28. Allen JD, Pernet B. Intermediate modes of larval development: bridging the gap between planktotrophy and lecithotrophy. Evol Dev. 2007;9:643–53.

    Article  PubMed  Google Scholar 

  29. Collin R. Nontraditional life-history choices: what can “intermediates” tell us about evolutionary transitions between modes of invertebrate development? Integr Comp Biol. 2012;52:128–37.

    Article  PubMed  Google Scholar 

  30. Zakas C, Rockman MV. Dimorphic development in Streblospio benedicti: genetic analysis of morphological differences between larval types. Int J Dev Biol. 2014;58:593–9.

    Article  PubMed  Google Scholar 

  31. Emlet RB. Facultative planktotrophy in the tropical echinoid Clypeaster rosaceus (Linnaeus) and a comparison with obligate planktotrophy in Clypeaster subdepressus (gray) (Clypeasterodia: Echinoidea). J Exp Mar Biol Ecol. 1986;95:183–202.

    Article  Google Scholar 

  32. Allen JD, Zakas C, Poldolsky RD. Effects of egg size reduction and larval feeding on juvenile quality for a species with facultative-feeding development. J Exp Mar Biol Ecol. 2006;331:186–97.

    Article  Google Scholar 

  33. Zigler KS, Lessios HA, Raff RA. Egg energetics, fertilization kinetics, and population structure in echinoids with facultatively feeding larvae. Biol Bull. 2008;215:191–9.

    Article  PubMed  Google Scholar 

  34. Armstrong AF, Lessios HA. The evolution of larval developmental mode: insights from hybrids between species with obligately and facultatively planktotrophic larvae. Evol Dev. 2015;17:278–88.

    Article  PubMed  Google Scholar 

  35. Reitzel AM, Miner BG. Reduced planktotrophy in larvae of Clypeaster rosaceus (Echinodermata, Echiniodea). Mar Biol. 2007;151:1525.

    Article  Google Scholar 

  36. Gildor T, Malik A, Sher N, Avraham L, de-Leon SBT. Quantitative developmental transcriptomes of the Mediterranean Sea urchin Paracentrotus lividus. Mar Genomics. 2016;25:89–94.

    Article  PubMed  Google Scholar 

  37. Gaitán-Espitia JD, Hofmann GE. Gene expression profiling during the embryo-to-larva transition in the giant red sea urchin Mesocentrotus franciscanus. Ecol Evol. 2017;7:2798–811.

    Article  PubMed  PubMed Central  Google Scholar 

  38. Wygoda JA, Yang Y, Byrne M, Wray GA. Transcriptomic analysis of the highly derived radial body plan of a sea urchin. Genome Biol Evol. 2014;6:964–73.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  39. Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, Adiconis X, Fan L, Raychowdhury R, Zeng Q, Chen Z, Mauceli E, Hacohen N, Gnirke A, Rhind N, di Palma F, Birren BW, Nusbaum C, Lindblad-Toh K, Friedman N, Regev A. Full-length transcriptome assembly from RNA-seq data without a reference genome. Nat Biotechnol. 2011;29:644–52.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  40. Li W, Godzik A. CD-hit: a fast program for clustering and comparing large sets of protein or nucleotide sequences. Bioinformatics. 2006;22:1658–9.

    Article  PubMed  CAS  Google Scholar 

  41. Haas BJ, Papanicolaou A, Yassour M, Grabherr M, Blood PD, Bowden J, Couger MB, Eccles D, Li B, Lieber M, Macmanes MD, Ott M, Orvis J, Pochet N, Strozzi F, Weeks N, Westerman R, William T, Dewey CN, Henschel R, Leduc RD, Friedman N, Regev A. De novo transcript sequence reconstruction from RNA-seq using the trinity platform for reference generation and analysis. Nat Protoc. 2013;8:1494–512.

    Article  PubMed  CAS  Google Scholar 

  42. Meyer E. Accessed 1 Feb 2017.

  43. Parra G, Bradnam K, Korf I. CEGMA: a pipeline to accurately annotate core genes in eukaryotic genomes. Bioinformatics. 2007;23:1061–7.

    Article  PubMed  CAS  Google Scholar 

  44. Simão FA, Waterhouse RM, Ioannidis P, Kriventseva EV, Zdobnov EM. BUSCO: assessing genome assembly and annotation completeness with single-copy orthologs. Bioinformatics. 2015;31:3210–2.

    Article  PubMed  CAS  Google Scholar 

  45. Emms DM, Kelly S. OrthoFinder: solving fundamental biases in whole genome comparisons dramatically improves orthogroup inference accuracy. Genome Biol. 2015;16:157.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  46. Cameron RA, Samanta M, Yuan A, He D, Davidson E. SpBase: the sea urchin genome database and web site. Nucleic Acids Res. 2008;37:D750–4.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  47. Conesa A, Götz S, García-Gómez JM, Terol J, Talón M, Robles M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005;21:3674–6.

    Article  PubMed  CAS  Google Scholar 

  48. Thomas PD, Kejariwal A, Campbell MJ, Mi H, Diemer K, Guo N, Ladunga I, Ulitsky-Lazareva B, Muruganujan A, Rabkin S, Vandergriff JA. PANTHER: a browsable database of gene products organized by biological function, using curated protein family and subfamily classification. Nucleic Acids Res. 2003;31:334–41.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  49. Ashburner M, Ball CA, Blake JA, Botstein D, Butler H, Cherry JM, Davis AP, Dolinski K, Dwight SS, Eppig JT, Harris MA. Gene ontology: tool for the unification of biology. Nat Genet. 2000;25:25.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  50. Falcon S, Gentleman R. Using GOstats to test gene lists for GO term association. Bioinformatics. 2006;23:257–8.

    Article  PubMed  CAS  Google Scholar 

  51. Robinson MD, Oshlack A. A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biol. 2010;11:R25.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  52. Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol. 2010;11:106.

    Article  CAS  Google Scholar 

  53. Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.

    Article  PubMed  CAS  Google Scholar 

  54. Gildor T, de-Leon SBT. Comparative study of regulatory circuits in two sea urchin species reveals tight control of timing and high conservation of expression dynamics. PLoS Genet 2015;11:e1005435.

  55. Samanta MP, Tongprasit W, Istrail S, Cameron RA, Tu Q, Davidson EH, Stolc V. The transcriptome of the sea urchin embryo. Science. 2006;314:960–2.

    Article  PubMed  CAS  Google Scholar 

  56. Wei Z, Angerer RC, Angerer LM. A database of mRNA expression patterns for the sea urchin embryo. Dev Biol. 2006;300:476–84.

    Article  PubMed  PubMed Central  CAS  Google Scholar 

  57. Tadros W, Lipshitz HD. The maternal-to-zygotic transition: a play in two acts. Development. 2009;136:3033–42.

    Article  PubMed  CAS  Google Scholar 

  58. Miner BG, McEdward LA, McEdward LR. The relationship between egg size and the duration of the facultative feeding period in marine invertebrate larvae. J Exp Mar Biol Ecol. 2005;321:135–44.

    Article  Google Scholar 

Download references


We are especially grateful to Haris Lessios for his feedback and guidance on this project. Additionally, we thank Rachel Collin, Plinio Gondola, and the rest of the staff at the Bocas del Toro research station for their help during our tissue collection in Panama. Brenda Cameron assisted with RNA extractions and library preparations. Serena Caplins provided bioinformatic consultation particularly with OrthoFinder.


AFA was supported by a short-term fellowship from the Smithsonian Tropical Research Institute, as well as graduate research fellowship from the National Science Foundation (DGE-1650042), and a research grant through the UC Davis Center for Population Biology. This research was also supported by a grant to RKG from the US National Science Foundation (OCE-1459815). The computational analyses for this project were run on the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number ACI-1548562.

Availability of data and materials

The final transcriptomes, blast2GO annotations, as well as the raw sequenced reads used for differential expression analysis, are available on DRYAD (

Author information

Authors and Affiliations



This project was conceived and designed by AFA and RKG. AFA conducted the experiments, RNA extractions, and RNAseq library preparations. AFA assembled the transcriptomes and conducted the data analyses presented here. RKG and AFA wrote and edited the manuscript together. Both authors read and approved the final manuscript.

Corresponding author

Correspondence to Anne Frances Armstrong.

Ethics declarations

Ethics approval and consent to participate

Collection and export permits for this project were granted to AFA from the Authority of the Aquatic Resources of Panama (ARAP permit NO 58, 59).

Competing interests

The authors declare they have no competing interests.

Publisher’s Note

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

Additional file

Additional file 1:

Table S1. Text table listing every gene ontology term identified as significantly over expressed for each species at each developmental stage along with their p-values. (PDF 52 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (, 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 ( applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Armstrong, A.F., Grosberg, R.K. The developmental transcriptomes of two sea biscuit species with differing larval types. BMC Genomics 19, 368 (2018).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: