For 40 years, paleontological studies of marine gastropods have suggested that species selection favors lineages with short-lived (lecithotrophic) larvae, which are less dispersive than long-lived (planktotrophic) larvae. Although lecithotrophs appeared to speciate more often and accumulate over time in some groups, lecithotrophy also increased extinction rates, and tests for state-dependent diversification were never performed. Molecular phylogenies of diverse groups instead suggested lecithotrophs accumulate without diversifying due to frequent, unidirectional character change. Although lecithotrophy has repeatedly originated in most phyla, no adult trait has been correlated with shifts in larval type. Thus, both the evolutionary origins of lecithotrophy and its consequences for patterns of species richness remain poorly understood. Here, we test hypothesized links between development mode and evolutionary rates using likelihood-based methods and a phylogeny of 202 species of gastropod molluscs in Sacoglossa, a clade of herbivorous sea slugs. Evolutionary quantitative genetics modeling and stochastic character mapping supported 27 origins of lecithotrophy. Tests for correlated evolution revealed lecithotrophy evolved more often in lineages investing in extra-embryonic yolk, the first adult trait associated with shifts in development mode across a group. However, contrary to predictions from paleontological studies, species selection actually favored planktotrophy; most extant lecithotrophs originated through recent character change, and did not subsequently diversify. Increased offspring provisioning in planktotrophs thus favored shifts to short-lived larvae, which led to short-lived lineages over macroevolutionary time scales. These findings challenge long-standing assumptions about the effects of alternative life histories in the sea. Species selection can explain the long-term persistence of planktotrophy, the ancestral state in most clades, despite frequent transitions to lecithotrophy.

Figure S1, 16S secondary structure model for Sacoglossa

Secondary structure model for a portion of the mitochondrial large ribosomal subunit rRNA (16S) gene for Elysia crispata; positions conserved across Mollusca are bolded.

Fig S1, 16S secondary structure model.pdf

Figure S2, Alignment of partial 16S gene sequences

Alignment of partial 16S gene sequences for all taxa, annotated according to secondary structure. Regions of ambiguous alignment (indicated by green shaded boxes marked “X”) were manually deleted prior to analyses.

Fig S2, full 16S alignment annotated.pdf

Figure S3, Annotated alignment of partial 28S rDNA sequences

Alignment of a portion of the nuclear large ribosomal subunit rRNA (28S) gene for all taxa, showing annotated secondary structure based on features conserved across Metazoa. Regions of ambiguous alignment (indicated by green shaded boxes marked “X”) were manually deleted prior to analyses.

Fig S3, 28S alignment annotated.pdf

Figure S4, Individual gene trees (a-d), BI consensus tree (e), and ML tree with two data partitions (f)

Supplementary phylogenetic trees, including individual ML gene trees for each locus based on partial sequences of the (a) mitochondrial COI (658 bp), (b) mitochondrial 16S rRNA (404 bp), (c) nuclear 28S rRNA (1392 bp), and (c) nuclear H3 (328 bp) genes; (e) the BI consensus phylogram; and (f) the ML phylogram and significant bootstrap support values when using two data partitions as per PartitionFinder.

Fig S4a-f, gene trees, BI consensus, ML alt partitioned.pdf

Figure S5, Collapsed phylogeny showing distribution of species richness and development modes

Distribution of species diversity among supported clades on the ML phylogeny of Sacoglossa. Traditional genera or genus-level clades were collapsed in proportion to known species richness, listed for each group. Species diversity and a pruned tree were used as inputs for ML tests of state-dependent speciation rates using the BiSSE model. Larval development mode for terminal taxa are coded as: P, planktotrophic; L, lecithotrophic; ?, unknown. The number of species with each development mode (or missing data) is given in parentheses for collapsed clades or terminals representing multiple taxa, in the same order (P,L,?). The sum of these values represents the total number of species known in a clade, including undescribed candidate species identified in the present work. Red branches indicate shifts in diversification detected by MEDUSA that are independent of larval development mode. Traditional families are named to the right of, and delimited by, colored horizontal bars; bolded family names are traditionally placed in superfamily Limapontioidea. Select higher clade names are given on branches.

Fig S5, BiSSE tree v4.pdf

Figure S6, Posterior distributions of rate parameters from dependent models of correlated trait evolution

Posterior distributions of parameters from dependent models of trait evolution, allowing ECY and development mode to evolve in a correlated manner. Transition rates are labeled as in Fig. 2. Z = percentage of time a rate was in the zero bin.

Fig S6, dep model rate dist.tif

Table S1, Sampled taxa and collection details

Species names, sample codes and collection details for sequenced taxa used in phylogenetic analyses. Blank cells reflect missing information for sequence data obtained from a public database.

Table_S1_sampled_taxa+collection_details_v6.doc

Table S2, Developmental character data for Sacoglossa

Developmental character data for Sacoglossa including larval development mode, presence/absence and pattern of extra-capsular yolk (ECY), mean egg diameter (± SD), and mean larval shell width (± SD) at hatching. Taxa are listed alphabetically by traditional family within higher clades, and then by binomial name within family. Generic names in quotation marks denote taxa that do not group with most other members of the genus to which they are traditionally assigned.

Table_S2_devel_data_v7_refs.xls

Table S3, NCBI accession numbers

Accession numbers from the National Center for Bioinformatics database for sequences generated and/or analyzed in this study.

Table_S3_accession_nos.xls

Table S4, BiSSE model fit using a pruned input tree and global estimate of unsampled taxa

Comparison of BiSSE models correcting for missing data with a pruned input tree and one overall estimate of the percentage of unsampled taxa (69%), with either (a) one rate of character change, or (b) rates of reversal (q10) constrained to be rare (<1%) relative to gains of lecithotrophy (q01).

Table S4, split BiSSE model fit, % missing taxa.doc

Table S5, BiSSE parameter estimates using a percentage of unsampled taxa to correct for missing data

Maximum-likelihood parameter estimates for BiSSE models, correcting for missing data using an overall percentage of unsampled taxa (69%). The ML tree (Fig. 3) was pruned of terminals missing character data, and used as the input tree for BiSSE analyses with three phylogenetic partitions across Sacoglossa, and either (A) one rate of character change, or (B) rates of reversal to planktotrophy (q10) constrained to be rare (<1%) relative to forward rates (q01). Estimates from the preferred model are bolded, with alternatives shown in descending order of AIC scores.

Table S5, split BiSSE params, % missing taxa.doc

Table S6, References cited in Table 4

References cited in Table 4, from which data on development modes were taken for select clades in Heterobranchia and Caenogastropoda.

Table S6, references v2.doc

Data from: Species selection favors dispersive life histories in sea slugs, but higher per-offspring investment drives shifts to short-lived larvae

Data files

Abstract

Figure S1, 16S secondary structure model for Sacoglossa

Figure S2, Alignment of partial 16S gene sequences

Figure S3, Annotated alignment of partial 28S rDNA sequences

Figure S4, Individual gene trees (a-d), BI consensus tree (e), and ML tree with two data partitions (f)

Figure S5, Collapsed phylogeny showing distribution of species richness and development modes

Figure S6, Posterior distributions of rate parameters from dependent models of correlated trait evolution

Table S1, Sampled taxa and collection details

Table S2, Developmental character data for Sacoglossa

Table S3, NCBI accession numbers

Table S4, BiSSE model fit using a pruned input tree and global estimate of unsampled taxa

Table S5, BiSSE parameter estimates using a percentage of unsampled taxa to correct for missing data

Table S6, References cited in Table 4

Data from: Species selection favors dispersive life histories in sea slugs, but higher per-offspring investment drives shifts to short-lived larvae

Data files

Abstract

Usage notes

Figure S1, 16S secondary structure model for Sacoglossa

Figure S2, Alignment of partial 16S gene sequences

Figure S3, Annotated alignment of partial 28S rDNA sequences

Figure S4, Individual gene trees (a-d), BI consensus tree (e), and ML tree with two data partitions (f)

Figure S5, Collapsed phylogeny showing distribution of species richness and development modes

Figure S6, Posterior distributions of rate parameters from dependent models of correlated trait evolution

Table S1, Sampled taxa and collection details

Table S2, Developmental character data for Sacoglossa

Table S3, NCBI accession numbers

Table S4, BiSSE model fit using a pruned input tree and global estimate of unsampled taxa

Table S5, BiSSE parameter estimates using a percentage of unsampled taxa to correct for missing data

Table S6, References cited in Table 4

Works referencing this dataset