Skip to main content

It’s not what it looks like: molecular data fails to substantiate morphological differences in two sea hares (Mollusca, Heterobranchia, Aplysiidae) from southern Brazil


Species of sea hares have been recognized traditionally based on morphological traits, mainly the radula, external coloration, and reproductive anatomy. However, recent studies have shown substantial color variation in some sea slug species. Molecular data have been successfully used to differentiate morphologically similar species of “opisthobranchs” and resolve questions on the taxonomic value of color. The objective of this paper is to use molecular data in an attempt to elucidate whether specimens of Aplysia brasiliana with distinct colorations and morphologies are actually the same species. To this end, DNA from 14 specimens of A. brasiliana was extracted, including five specimens identified as a distinct morphotype from typical A. brasiliana. Although the two morphotypes have consistent differences in their external morphology and radula, the molecular data confirmed that there are no significant genetic differences between them. This is another example of the need to re-evaluate taxonomic decisions based on morphology in light of molecular evidence.


Species of sea hares have been recognized traditionally based on morphological traits, mainly the radula, external coloration, and reproductive anatomy (Eales 1960). The original descriptions of most species of Aplysia were based solely on external morphological characteristics, body and shell shape and particularly color pattern (see Eales 1960). Marcus and Marcus (1958) and Eales (1960) were the first authors to described species of Aplysia including illustrations of the radula and the reproductive anatomy, but always stressing color differences between species for taxonomic decisions.

The wide range of variation in the coloration of some species of Aplysia has been the source of considerable confusion as how to differentiate species (Winkler 1959). MacMunn (1899) and Winkler (1959) have shown that in Aplysia the normal skin pigment is basically made up of degradation products of the tetrapyrrole molecule of chlorophyll, in other words, there should be a close connection between diet and coloration. However, new species names for Aplysia and other “opisthobranch” molluscs have been introduced based on color differences with closely related species. Recent studies have shown substantial color variation in some “opisthobranch” species (Valdés et al. 2013) including sea hares (Alexander and Valdés 2013). The implication of these findings is that in Aplysia color differences among individuals may have limited taxonomic value.

In recent years, molecular data have been successfully used to differentiate morphologically similar species of “opisthobranch” sea slugs and to resolve taxonomic questions on the significance of color variation for systematics (Medina et al. 2001, 2005; Ornelas-Gatdula et al. 2011; Alexander and Valdés 2013). In this paper, we study two different morphotypes of Aplysia brasiliana (Rang, 1828) collected from Southern Brazil. These two color forms display consistent morphological differences, but the two morphotypes are found sympatrically, in the same habitat during the same time of the year. We use molecular data in an attempt to elucidate whether specimens with distinct colorations and morphologies are actually the same species.

Materials and methods

Source of material

DNA from 14 specimens of A. brasiliana was extracted (Table 1). Five specimens were identified as a distinct morphotype from typical A. brasiliana, here referred to as Aplysia aff. brasiliana. The specimens were collected, from November 2011 to January 2013, in Ponta da Praia, Santos, São Paulo, Brazil (23°59′29″S; 46°18′15″W) and in São Sebastião, São Paulo, Brazil (23°49′43.06″S; 45°25′21.83″W). Additionally, sequences of A. brasiliana from Florida, Aplysia parvula Mörch, 1863 from the Bahamas, Aplysia kurodai (Baba, 1937) from Japan, Aplysia dactylomela (Rang 1828) from the Bahamas and Aplysia argus Rüppell and Leuckart 1828 from the Hawaiian Islands ware included in this study for comparison using sequences available in GenBank. Phyllaplysia taylori Dall 1900 was used as the outgroup (Table 1).

Table 1 List of specimens included in this study Including Locality, Voucher Numbers and GenBank Accession Numbers

All the specimens from Brazil are deposited at the Museu de Zoologia Universidade de São Paulo (MZUSP). Additional material was obtained from the collections of the Department of Biological Sciences of the California State Polytechnic University, Pomona (CPIC) and Museu Nacional Universidade Federal do Rio de Janeiro (MNRJ).

Morphological study

All the specimens of A. brasiliana and A. aff. brasiliana were photographed alive in the field using a Nikon P100 digital camera. Several traits were recorded photographically, including the shell, the opening of the siphon and the gill, and the anterior end of the animal. All sea hares were photographed with the parapodia opened and closed. The specimens were taken to the laboratory at MZSP. All the individuals were euthanized by asphyxiation (kept in water without oxygenation) and preserved in 70 % EtOH. Under a stereomicroscope, the buccal mass was dissected from two specimens of each morphotype (A. brasiliana: MZSP 109978, MZSP 103251 and A. aff. brasiliana: MZSP 114920, MNRJ 18228). Tissue surrounding the radula and jaws was dissolved in 10 % sodium hydroxide for 2 days. The radulae were rinsed in water, dried, mounted on stubs, and coated with a gold–palladium alloy for scanning electron microscopy. The samples were observed under a scanning electron microscope (SEM) Zeiss LEO-440 at the Laboratório de Microscopia Eletrônica do Museu de Zoologia, Universidade de São Paulo.

DNA extraction

DNA was extracted from a small foot sample using a hot Chelex® protocol with modifications. The tissue was rinsed and rehydrated using 1.0 mL TE buffer (10 mM Tris, 1 mM EDTA, pH 8.0) for 20 min. A 10 % (w/v) Chelex® 100 (100–200 mesh, sodium form; Bio-Rad) was prepared using TE buffer. After rehydration, the mixture was then centrifuged, 975 mL of the supernatant was removed, and 175 mL of the Chelex® solution was added. Samples were then heated in a 56 °C water bath for 3 h, then heated in a 100 °C heating block for 8 min, and the supernatant was used for PCR.


Palumbi’s universal 16S rRNA primers (16S ar-L 5′-CGCCTGTTTATCAAAAACAT-3′, 16S br-H 5′-CCGGTCTGAACTCAGATCACGT-3′) (Palumbi 1996), Folmer’s universal cytochrome oxidase subunit 1 primers (LCO1490 5′-GGTCAACAAATCATAAAGATATTGG-3′, HCO2198 5′-TAAACTTCAGGGTGACCAAAAAATCA-3′) (Folmer et al. 1994) and Colgan’s universal Histone-3 primers (H3 AF 5′-ATGGCTCGTACCAAGCAGACGGC-3′, H3 AR 5′-ATATCCTTGGGCATGATGGTGAC-3′ developed by Colgan et al. 1998) were used to amplify the regions of interest for all specimens.

The internal primers 16Sar-FAP 5′-AAAGACGAGAAGACCCTTAGAGTTTT-3′ and 16Sbr-FAP 5′-AAAACTCTAAGGGTCTTCTCGTCTTT-3′ (Ornelas-Gatdula et al. 2011) were used to amplify 16S in the São Sebastião specimens (MZSP103233; MZSP103234; MZSP103250) and only partial sequences were obtained.

PCR amplification and sequencing

The master mix was prepared using 34.75 mL H2O, 5 mL PCR Buffer (500 mM KCl; 100 mM Tris, pH 8.3), 5 mL 50 mM MgCl2, 1 mL 40 mM dNTPs, 1 mL 10 mM primer 1, 1 mL 10 mM primer 2, 0.25 mL 5 mg/mL Taq, and 2 mL extracted DNA. Reaction conditions for H3 and 16S rRNA were as follows: an initial denaturation for 2 min at 94 °C, 35 cycles for 30 s at 94 °C, annealing for 30 s at 50 °C, and elongation for 1 min at 72 °C, and a final elongation for 7 min at 72 °C. Reaction conditions for CO1 (universal) were as follows: lid heated to 105 °C and initial denaturation for 4 min at 95 °C, 40 cycles for 45 s at 94 °C, annealing for 45 s at 48 °C, and elongation for 90 s at 72 °C, followed by a final elongation step for 10 min at 72 °C.

Gel electrophoresis was used to determine whether PCR resulted in products with the band size of the approximate length for the targeted gene (approximately 375 base pairs [bp] for H3; 475 bp for 16S rRNA, 195 bp for 16S fragments, 700 bp for CO1), then were purified using the GeneJET PCR Purification Kit (Thermo Scientific). Cleaned PCR samples were quantified using a NanoDrop 1000 Spectrophotometer (Thermo Scientific). These samples were sequenced at the Eton Biosciences (San Diego, CA).

Sequence analyses

Sequences for each gene were assembled and edited using Geneious Pro 4.7.4 (Drummond et al. 2009). Geneious was also used to extract the consensus sequence and to construct the alignment for each gene using the default parameters and to concatenate the alignments. The sequences were not trimmed after alignment. The lengths of the gene fragments used in the phylogenetic analyses are 305–328 bp for H3, 440–454 bp for 16S, and 501–658 bp for COI. To assess whether the three genes have significantly conflicting signals, the incongruence length difference (ILD) test (Mickevich and Farris 1981; Farris et al. 1994), implemented in PAUP*4.0 as the partition homogeneity test (Swofford 2002), was calculated for all genes combined with 2,000 replicates. Phylogenetic analyses were conducted for all genes concatenated including those specimens for which at least two genes (including H3) were available. The best-fit models of evolution (GTR + I for CO1; HKY + I for 16S, and GTR + I+G for H3) were determined using the Akaike information criterion (Akaike 1974) implemented in MrModelTest (Nylander 2004). Bayesian analyses were conducted using MrBayes 3.2 (Ronquist et al. 2012), partitioned by gene (unlinked). The Markov chain Monte Carlo analysis was run with two runs of six chains for 20 million generations, with sampling every 100 generations. The default 25 % burn-in was applied before constructing the majority-rule consensus tree. Convergence was inspected using Tracer 1.5 (Rambaut and Drummond 2007). Maximum-likelihood analysis was conducted for the entire concatenated alignment with RAxML GUI 1.0 (Silvestro and Michalak 2012) using the bootstrap + consensus option (10,000 replicates) and the GTR + I model.

To visualize the genetic structure of A. brasilana, a haplotype network was constructed using all COI sequences of A. brasiliana using the program TCS 1.21 (Clement et al. 2000) with a 95 % connection limit. The same analyses were not conducted for 16S and H3 because all the sequences for A. brasiliana and A. aff brasiliana were identical.

Arlequin 3.5 (Excoffier and Lischer 2010) was also used to calculate FST values as a measure of pairwise differences between both morphotypes of A. brasiliana for both COI and 16S. The significance of the pairwise FST value was estimated by performing 10,000 permutations.


Morphological study

The morphological examination of the specimens of A. brasiliana from Brazil revealed the presence of two distinct morphotypes. The main differences found between morphotypes include the external coloration and morphology. Specimens of morphotype 1 (the typical Aplysia brasiliana) posses a greenish and/or brownish background color, with few yellowish patches on the body and parapodia (Fig. 1a); the parapodia are relatively narrow, as wide as the rest of the body; there is a small pore (p) above the shell, which is covered by a papilla (Fig. 1b); the edge of the siphon does not completely cover the opening (s1), forming a semicircle (Fig. 1c); the oral tentacles are small (Fig. 1d). Morphotype 2 (Aplysia aff. brasiliana) has a black, gray or dark reddish background color, with conspicuous white spots scattered throughout the mantle and parapodia (Fig. 1e); the parapodia are broad, twice as wide as the rest of the body; the small pore above the shell is not covered by a papilla (Fig. 1f); the edge of the siphon covers the opening completely (s2) forming a circle (Fig. 1g); and the oral tentacles are large (Fig. 1h).

Fig. 1

a–d Aplysia brasiliana morphotype. a Background greenish and/or brownish with few yellowish patches around the body and parapodium. b Parapodium end a small pore (p) covered by a papilla. c The edge of the siphon (s1). d Oral tentacles. e–h Aplysia aff brasiliana morphotype. e Background black, gray or dark reddish with predominant white spots scattered throughout the mantle and parapodium. f Parapodia. g The edge of the siphon (s2). h Oral tentacles

Examination of the radular morphology also revealed some differences between A. brasiliana and A. aff brasiliana. The rachidian teeth of A. brasiliana are broad with a denticulate central cusp and a secondary smaller cusp on each side; the bases of the teeth are convex (Fig. 2a). The mid-lateral teeth have an elongate cusp with several denticles on each side (Fig. 2b). The outer teeth are similar in shape (Fig. 2c). On the contrary, A. aff. brasiliana has much narrower rachidian teeth, with a shorter, more rounded central cusp, the teeth are triangular in shape, and the bases are less convex than in A. brasiliana (Fig. 2d). The lateral teeth have shorter cusps than in A. brasiliana (Fig. 2e). The outermost lateral teeth are very similar to those of A.brasiliana (Fig. 2f).

Fig. 2

a–c Aplysia brasiliana radula. a Rachidian tooth. b Lateral teeth. c Marginal teeth. d–f A. aff brasiliana radula. d Rachidian tooth. e Lateral teeth. f Marginal teeth

Molecular study

The ILD test showed no significant conflicting signals between the three genes combined 16S versus CO1 (P = 0.918), CO1 versus H3 (P = 0.909), and 16S versus H3 (P = 1). The maximum-likelihood bootstrap consensus and Bayesian consensus trees for all three genes combined recovered three clades (Fig. 3). One clade, including all specimens of A. brasiliana and A. aff brasiliana had a posterior probability (PP) of 0.88 and bootstrap value (BS) of 100. A second clade including specimens of A. dactylomela PP of 1 and BS of 98. A third clade including specimens of A. kurodai from Japan had a PP of 1 and a BS of 100. The specimens of A. brasiliana and A. aff brasiliana do not form distinct clades.

Fig. 3

Schematic representation of the Bayesian consensus phylogenetic tree of Aplysiidae based on three molecular markers. Statistical support percentages are shown on internal branches. Bayesian posterior probabilities precede maximum-likelihood bootstrap values

The CO1 haplotype network (Fig. 4) reveals the presence of 6 distinct haplotypes in A. brasiliana including the specimen from Florida. However, there is no clear structure in the haplotype network in relation to the morphology of the animals, and one of the sequences of A. aff. brasiliana is identical to three sequences of A. brasiliana.

Fig. 4

Haplotype network (with sequence IDs) of CO1 haplotypes generated with TCS. Pie charts indicate the proportion of different haplotype clusters found in each locality. Haplotypes identified from the two morphotypes of A. brasiliana from Brazil and from other localities are color coded

The result of the pairwise FST analysis between the two morphotypes was not significant for both COI and 16S (Table 2), confirming that there are no significant genetic differences between the two morphotypes.

Table 2 Population comparisons results for COI and 16S, with FST values and associated p values


Although color pattern and morphological traits have been used traditionally to distinguish “opisthobranch” sea slug species, recent studies have challenged this view in light of newly available molecular data (Pola et al. 2006; Ornelas-Gatdula et al. 2011; Alexander and Valdés 2013; Valdés et al. 2013). The present study shows that specimens of A. brasiliana with distinctive morphologies and colorations are not genetically different and constitute morphotypes of the same species. As in some of the other examples mentioned above, the differences between the two morphotypes of A. brasiliana appear to be discrete and consistent, with no apparent intermediate specimens between the two morphotypes and detectable in a range of internal and external traits.

The radular differences are particularly conspicuous. Other studies on sea hares have shown that some species display substantial radular variability. For example, Martínez and Ortea (2002) showed that radular morphology of Aplysia fasciata can vary between large and small animals. Although this has not been tested for sea hares, it has been shown that diet can influence tooth morphology (size and number) in other “opisthobranch” sea slugs. For example, in sacoglossans, smooth teeth can be associated with filamentous food (Bleakney 1990; Jensen 1993). Bleakney (1989) also noted some other morphological differences in the teeth of a species of “opisthobranchs” from different geographic regions.

The combination of morphological and molecular data in an integrative framework is the ideal approach to resolve difficult taxonomic problems. This is particularly important in “opisthobranchs”, which appear to be morphologically plastic. Because of this and other studies, it is becoming more evident that many taxonomic decisions based on morphology will have to be re-evaluated in light of molecular evidence. Radular morphology and external coloration have been widely used to distinguish species of “opisthobranch”, but these data may not be as useful as previously thought. This hypothesis will have to be tested in a wider variety of taxa.


  1. Akaike H (1974) A new look at the statistical model identifications. IEEE Trans Autom Control 19:716–723

    Article  Google Scholar 

  2. Alexander J, Valdés A (2013) The ring doesn’t mean a thing: molecular data suggest a new taxonomy for two Pacific species of sea hares (Mollusca: Opisthobranchia, Aplysiidae). Pac Sci 67:283–294. doi:10.2984/67.2.10

    Article  Google Scholar 

  3. Bleakney JS (1989) Morphological variation in the radula of Placida dendritica (Alder and Hanrork, 1843) (Opisthobranchia: Ascoglossa/Sacoglossa) from Atlantic and Pacific populations. Veliger 32:171–181

    Google Scholar 

  4. Bleakney JS (1990) Indirect evidence of a morphological response in the radula of Placida dendritica (Alder and Hancock, 1843) (Opisthobranchia: Ascoglossa/Sacoglossa) to different algal prey. Veliger 33:111–115

    Google Scholar 

  5. Clement M, Posada D, Crandall KA (2000) TCS: a computer program to estimate gene genealogies. Mol Ecol 9:1657–1659. doi:10.1046/j.1365-294x.2000.01020.x

    CAS  PubMed  Article  Google Scholar 

  6. Colgan DJ, McLauchlan A, Wilson GDF, Livingston SP, Edgecombe GD, Macaranas J, Cassis G, Gray MR (1998) Histone H3 and U2 snRNA DNA sequences and arthropod molecular evolution. Aust J Zool 46:419–437. doi:10.1071/ZO98048

    Article  Google Scholar 

  7. Drummond AJ, Ashton B, Cheung M, Heled J, Kearse M, Moir R, Stones-Havas S, Thierer T, Wilson AC (2009) Geneious ver. 4.8.

  8. Eales NB (1960) Revision of the world species of Aplysia (Gastropoda Opisthobranchia). Bull Brit Mus (Nat Hist) Zool 5:267–404

  9. Excoffier L, Lischer HEL (2010) Arlequin suite ver 3.5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol Ecol Resour 10:564–567. doi:10.1111/j.1755-0998.2010.02847.x

    PubMed  Article  Google Scholar 

  10. Farris JS, Kallersjo M, Kluge AG, Bult C (1994) Testing significance of incongruence. Cladistics 10:315–319. doi:10.1111/j.1096-0031.1994.tb00181.x

    Article  Google Scholar 

  11. Folmer O, Black M, Hoeh W, Lutz R, Vrijenhoek R (1994) DNA primers for amplification of mitochondrial cytochrome c oxidase subunit I from diverse metazoan invertebrates. Mol Mar Biol Biotechnol 3:294–299

    CAS  PubMed  Google Scholar 

  12. Jensen KR (1993) Morphological adaptations and plasticity of radular teeth of the Sacoglossa (= Ascoglossa) (Mollusca: Opisthobranchia) in relation to their food plants. Biol J Linn Soc 48:135–155

  13. MacMunn CA (1899) The pigments of Aplysia punctata. J Physiol 24:1–10

    CAS  Google Scholar 

  14. Marcus E, Marcus E (1958) Notes on Aplysia. Univ São Paulo, Bull Inst Oceanogr 8:3–22

    Article  Google Scholar 

  15. Martínez E, Ortea J (2002) On the synonymy between Aplysia winneba Eales, 1957 and Aplysia fasciata Poiret 1789 (Mollusca: Opisthobranchia: Anaspidea). Iberus 20:11–21

    Google Scholar 

  16. Medina M, CollinsTM Walsh PJ (2005) Phylogeny of sea hares in the Aplysia clade based on mitochondrial DNA sequence data. Bull Mar Sci 76:691–698

    Google Scholar 

  17. Medina M, Collins TM, Walsh PJ (2001) mtDNA ribosomal gene phylogeny of sea hares in the genus Aplysia (Gastropoda, Opisthobranchia, Anaspidea): implications for comparative neurobiology. Syst Biol 50:676–688. doi:10.1080/106351501753328802

    CAS  PubMed  Article  Google Scholar 

  18. Mickevich MF, Farris JS (1981) The implications of congruence in Menidia. Syst Zool 30:351–370

    Article  Google Scholar 

  19. Nylander JAA (2004) MrModeltest ver.2. Evolutionary Biology Centre, Uppsala University, Uppsala.

  20. Ornelas-Gatdula E, Dupont A, Valdés A (2011) The tail tells the tale: taxonomy and biogeography of some Atlantic Chelidonura (Gastropoda: Cephalaspidea: Aglajidae) inferred from nuclear and mitochondrial gene data. Zool J Linn Soc 163:1077–1095. doi:10.1111/j.1096-3642.2011.00749.x

    Article  Google Scholar 

  21. Palumbi SR (1996) Nucleic acids II: the polymerase chain reaction. In: Hillis DM, Moritz C, Mable BK (eds) Molecular systematics. Sinauer, Sunderland, pp 205–247

    Google Scholar 

  22. Pola M, Vallès Y, Cervera JL, Medina M, Gosliner TM (2006) Taxonomic status of Tambja abdere and Tambja fusca based on morphological and molecular evidence, with comments on the phylogeny of the subfamily Nembrothinae (Nudibranchia: Polyceridae). Ann Zool Fenn 43:52–64

    Google Scholar 

  23. Rambaut A, Drummond AJ (2007) Tracer v1.5.

  24. Ronquist F, Teslenko M, Van der Mark P, Ayres DL, Darling A, Hohna S, Larget B, Liu L, Suchard MA, Huelsenbeck J (2012) MrBayes 3.2: Efficient Bayesian phylogenetic inference and model choice across a large model space.

  25. Silvestro D, Michalak I (2012) Raxml GUI: a graphical front-end for RAxML. Org Divers Evol 12:335–337. doi:10.1007/s13127-011-0056-0

    Article  Google Scholar 

  26. Swofford DL (2002) PAUP*: Phylogenetic analysis using parsimony (*and other methods). Ver. 4. Sinauer, Sunderland, Massachusetts

  27. Valdés A, Ornelas-Gatdula E, Dupont A (2013) Color pattern variation in a shallow-water species of Opisthobranch Mollusc. Biol Bull 224:35–46

    PubMed  Google Scholar 

  28. Winkler LR (1959) A mechanism of color variation operating in the west coast sea hare, Aplysia californica Cooper. Pac Sci 13:63–66

    Google Scholar 

Download references


We are grateful to Dr. Alexandre Pimenta (MNRJ) for loan of specimens for this project and to Dr. Daniela Sudatti (Universidade Federal Fluminense, UFF) for the donation of part of material here studied. We are grateful to Lara Guimarães, Laboratório de Microscopia Eletrônica, MZSP, for assisting with SEM examination of the specimens. Jermaine Mahuib helped to train LS in molecular work techniques. This project was supported by FAPESP (Fundação de Amparo à Pesquisa do Estado de São Paulo), proc. #2013/01000-4 and #2011/01493-5 to the senior author, procs. #2010/11253-9 to second author and a Cal Poly Pomona Provost Teacher–Scholar Award to the junior author.

Author information



Corresponding author

Correspondence to Luiza de Oliveira Saad.

Additional information

Communicated by H.-D. Franke.

Rights and permissions

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Saad, L.O., Cunha, C.M., Colpo, K.D. et al. It’s not what it looks like: molecular data fails to substantiate morphological differences in two sea hares (Mollusca, Heterobranchia, Aplysiidae) from southern Brazil. Helgol Mar Res 68, 523–530 (2014).

Download citation


  • Aplysia
  • CO1
  • 16S
  • Morphological variation