Phylogenomics illuminates the evolution of bobtail and bottletail squid (order Sepiolida)

Genome skimming provides robust phylogeny

Pioneering molecular phylogenetic studies in Sepiolida that used short regions of a few mitochondrial and nuclear genes failed to resolve the relationship of major clades^9,22,23. To increase the number of phylogenetically informative sites, Sanchez et al.¹¹ sequenced and analyzed the transcriptomes of multiple species of Euprymna Steenstrup, 1887, related bobtail squids including Sepiola parva Sasaki, 1913 and Sepiola birostrata Sasaki, 1918, and several bottletail squids. They found that S. parva grouped with the Euprymna species to the exclusion of S. birostrata, and further morphological analysis led to the formal redefinition of the genus Euprymna and the reassignment of S. parva Sasaki, 1913 to Euprymna parva¹¹. The following year, in an exhaustive study of hectocotylus structure, Bello¹⁹ proposed that Euprymna be split back into the original Euprymna Steenstrup, 1998 and a newly defined genus, Eumandya Bello 2020 that contains E. parva Sasaki, 1913 and E. pardalota Reid 2011, two taxa whose arms have two rows of suckers rather than four as in other Euprymna species. Similarly, Bello introduced a new genus, Lusepiola Bello, 2020 that has the effect of renaming Sepiola birostrata Sasaki, 1918 as Lusepiola birostrata. For clarity, we adopt the finer-grained nomenclature of Bello below, but happily note that E. parva and E. pardalota have the same abbreviations in both the notation of Sanchez et al.¹¹ and Bello¹⁹.

Sanchez et al.¹¹ also emphasized the need for more taxon sampling, careful species assignment, and the inclusion of more informative sites when studying this group of cephalopods. However, the distribution and lifestyle of many lineages of Sepiolida makes the collection of fresh tissue for RNA sequencing very challenging. To overcome this limitation, we sequenced the genomic DNA of several Sepiolida species at shallow coverage up to 3.6× and accessed by this way several mitochondrial and nuclear loci. Most of our samples were carefully identified at the species level based on morphological characters.

We recovered the mitochondrial genomes of the species targeted in this study and annotated the 13 protein-coding genes, 22 tRNAs, and two rRNAs (although only the conserved region of the large and small rRNA was obtained for Rondeletiola minor Naef, 1912).

Additionally, we also downloaded the complete mitochondrial genomes of S. austrinum and Idiosepius sp., and the transcriptome of E. tasmanica available in the NCBI database. The transcriptome of E. tasmanica was used to extract its complete set of mitochondrial protein-coding genes. We could reconstruct the mitochondrial gene order for all species with complete mtDNA genomes, but we observed no re-arrangement for members of Sepiolidae, and only Sepiadarium austrinum deviated from the arrangement seen in all other Sepiadariidae (Fig. S1).

To complement the mitochondrial-based evolutionary history, we also annotated several nuclear loci. As ribosomal gene clusters are present in numerous copies, they were successfully retrieved for almost all the species, except for 28S of the Sepiadariidae sp. specimen, which appeared problematic and was excluded.

By mapping reads to the reference genome of E. scolopes, we obtained 3,279,410 loci shared between at least two species and further selected 5215 loci presented in most of our Sepiolidae species, but allowed some missing data in the Euprymna + Eumandya clade. This was done because the phylogenetic relationships of the Euprymna + Eumandya species were previously described in detail in Sanchez et al.¹¹ using transcriptome data. Out of the 5215 loci, 5164 loci had a per-site coverage ranging between two and five. After trimming and removal of regions without informative sites, 577 loci remained. These ultraconserved loci had lengths ranging between 10 and 690 base pairs (bp), with an average of 65 bp. Our alignment matrix had a length of 37,512 bp and consisted of 16,495 distinct site patterns, and variable sites between 1 and 130 bp with an average value of 7 bp. We expected a low value of variable sites because these regions are highly conserved.

We considered resolved nodes to be those with the ultrafast bootstrap support and posterior probability larger than 95% and 0.9, respectively. Only the very unresolved nodes were found based on the mito_nc matrix (Fig. 1). However, among the species in these nodes, Adinaefiola ligulata Naef, 1912 was well supported with amino acid sequences from mitochondrial genes (posterior probability of 1 and 94% bootstrap support) and partially by the ultraconserved loci (posterior probability of 1, but only 85% bootstrap support) as sister to the Sepiola clade (Figs. 2 and S2). Moreover, compared to the mito_nc matrix and with identical topology, mito_aa and UCEbob fully resolved the relationship of the Indo-Pacific and Mediterranean Sea Sepiolinae. The tree generated by the nuclear_rRNA produces a topology with most nodes unsupported (Fig. S3), suggesting these markers are too conserved for assessing the relationships among this group.

Using the UCEbob matrix, the topology and supported relationships of Euprymna + Eumandya species resemble those reported in Sanchez et al.¹¹ using transcriptome sequences, proving our protocol valid when using low coverage sequencing and when a reference genome of the closest related species is available.

The position of R. minor showed discordance between mitochondrial and nuclear datasets. Using the mitochondrial matrices, R. minor rendered the Sepietta Naef, 1912 clade paraphyletic, whereas using the UCEbob and rRNA_nc matrices, R. minor appeared sister to the Sepietta clade. These relationships were resolved in both mitochondrial and nuclear-based trees and require further investigation with more DNA markers and a wider population sampling.

Molecular systematics of Sepiolida clades

Using the complete mitochondrial genome, ribosomal nuclear genes, and ultraconserved loci, we recovered the monophyly of the two families of the order Sepiolida—Sepiadariidae and Sepiolidae^9,24—and the monophyly of the three described subfamilies of the family Sepiolinae. However, contrary to what is proposed based on morphology in Young²⁴, the Rossinae is not sister to all the remaining Sepiolidae but rather is sister to Heteroteuthinae, although this is unresolved in the UCE phylogeny. With the lack of systematic work on these subfamilies, our robust phylogenetic backbone in Sepiolida using new samples carefully identified by morphology and with museum vouchers, represents a notable advance to clarify the evolution of morphological traits in major clades within the family.

Based on morphological characters of the hectocotylus, Bello¹⁹ recently split the polyphyletic Sepiola Leach 1817 into Lusepiola, Adinaefiola, and Boletzkyola, reserving Sepiola for the S. atlantica group sensu Naef 1923. These newly defined clades are consistent with our molecular phylogeny here and in Sanchez¹¹, who also noted the polyphyly of Sepiola in the Indo-Pacific lineage.

We find that Sepiolinae can be robustly split into two geographically distinct tribes: one that comprises species with known distribution in the Indo-Pacific region (tribe Euprymmini new tribe, defined as Sepiolinae with a closed bursa copulatrix, type genus Euprymna) and the other including all the Mediterranean and Atlantic species (tribe Sepiolini Appellof, 1989, defined here as Sepiolinae with an open bursa copulatrix, type genus Sepiola). Our molecular relationship is consistent with 13 of the 15 apomorphies used in the cladogram shown in Fig. 21 in Bello¹⁹. The other two proposed apomorphies in Bello (his apomorphic characters 4 and 6) group two IP lineages, Lusepiola and Inioteuthis, in a clade with species from the Mediterranean and Atlantic. Such relationships contradict our Euprymmini-Sepiolini sister relationship. Moreover, according to our phylogeny, apomorphy 6 of Bello, characterized by the participation of ventral and dorsal pedicels in the formation of the hectocotylus copulatory apparatus, implies that the male ancestor of Sepiolinae had a more developed hectocotylus that was simplified in the Euprymna and Eumandya clades.

Among euprymins, we confirmed the monophyly of Euprymna Steenstrup 1887 as found previously by transcriptome analysis¹¹. We also support the monophyly of Eumandya Bello, 2020 (Figs. 1 and 2), grouping the type species E. pardalota with E. parva along with the unnamed “Type 1” Ryukyuan species of Sanchez et al.¹¹, for which only hatchlings were available. The phylogenomic grouping of Ryukyuan “Type 1” with Eumandya suggests that when its adults are found (or hatchlings are raised to maturity), its arms will carry two rows of suckers. We also found an adult of a Ryukyuan “Type 4” (extending the notation of Sanchez et al.¹¹ in the coastal waters of Kume Island, that groups with E. scolopes from Hawaii, suggesting a divergence based on geographic isolation in the North Pacific. We also find that Lusepiola birostrata (formerly Sepiola birostrata) is grouped with Inioteuthis japonica as sister to a clade containing Euprymna, Eumandya, and an unnamed sepioline from Port Kembla, at the northeast of Martin Island in Australian waters.

Among the sepiolins, we confirm the monophyly of Sepietta (only for nuclear-genome-based trees, see below). Adinaefiola, another genus erected by Bello¹⁹, with Sepiola ligulata Naef 1912 as its type species; was found sister to the Sepiola clade, but only in the tree based on amino acid mitochondrial sequences (mito_aa matrix) with a bootstrap value of 94% and a posterior probability of 1 (Fig. S2).

Outside the sepiolines, members of the subfamily Heteroteuthinae are the most elusive and underrepresented in studies of cephalopod systematics due to their oceanic lifestyles. The placement of several heteroteuthin remains controversial. Lindgren et al.⁹, with six nuclear and four mitochondrial genes downloaded from GenBank found that Sepiolina Naef, 1912 was sister to Heteroteuthis Gray, 1849 + Rossia Owen, 1834+ Stoloteuthis Verril, 1881; rendering the subfamily Heteroteuthinae polyphyletic. In contrast, our work supports the monophyly of Heteroteuthinae by including Stoloteuthis and Heteroteuthis in this subfamily, while Rossia was placed within the Rossinae (Figs. 1, 2, S2). Members of Heteroteuthinae included in this study formed a sister group to a monophyletic Rossinae (Figs. 1, 2, S2). Semirossia, however, rendered the Rossia clade paraphyletic. Further discussion about the position of Semirossia is difficult because of the lack of information about the original source of this specimen in Kawashima et al.²⁵.

The light organ and luminescence evolution

Bobtail squids are thought to use the bioluminescence of their light organ to camouflage them from predators while foraging and swimming at night through a mechanism called counter-illumination. This has been researched extensively using E. scolopes as a model system^26,27,28. Unfortunately, the limited number of sequences available and the misidentification of bobtail squids in the GenBank database^11,29,30 have hindered our understanding of the light organ evolution in the whole taxon.

Our robust phylogeny and Bayesian reconstruction of ancestral bioluminescence clarify how the light organ and its luminescence have evolved in the family Sepiolidae. Members of Sepiolinae comprise neritic and benthic adults with bilobed light organs, except for two genera: Inioteuthis from the Indo-Pacific region, and the Sepietta species from the Mediterranean Sea and the Atlantic waters. The ancestor of the Sepiolinae very likely possessed a bilobed light organ that harbored luminescent symbiotic bacteria (Fig. 3). This character persisted until the ancestor of the euprymnins and sepiolins. Assuming that R. minor is sister to the Sepietta clade (as shown with the nuclear-based dataset, Fig. 2), it is clear that the bilobed light organ was lost once in Inioteuthis and Sepietta, and simplified to a rounded organ in R. minor. The alternative scenario, where R. minor renders the Sepietta clade paraphyletic (based on mitochondrial matrices, Fig. 1), is less plausible as it implies that the light organ was lost twice in the Sepietta group, once in S. obscura and then in the ancestor of S. neglecta and S. oweniana; or alternatively that it was lost in the ancestor of Sepietta–Rondelentiola followed by a reversion of this character in the lineage of Rondelentiola.

The light organ is also present in all members of Heteroteuthinae. These bobtails are pelagic as adults, and their light organ appears as a single visceral organ rather than the bilobed form found in nektobenthic Sepiolinae. In contrast to the bacteriogenic luminescence of the light organ in E. scolopes³¹, previous studies in H. dispar³ failed to detect symbiotic bacteria and suggested that the luminescence has an autogenic origin. Thus, it seems plausible that the monophyly of Heteroteuthinae found in our study supports the findings in Lindgren et al.⁹ for convergent evolution of autogenic light organs associated with pelagic lifestyle in many squid, octopus, and Vampyroteuthis Chun, 1903^9,32.

Divergence time of Sepiolida

The absence of fossils for this group limited our calculations of divergence time to the use of secondary calibrations. These calibrations can provide more accurate estimates depending on the type of primary calibrations that are used³³. We retrieved secondary calibrations from previous estimations in Tanner et al.¹⁵, who used eleven fossil records spanning from coleoids to gastropods in transcriptome-based phylogenetic trees. Specifically, we used the time for the splits of Sepia esculenta and S. officinalis (~91 Mya), Idiosepiidae, and Sepiolida (~132 Mya) and the origin of the Decapodiformes (root age, ~174 Mya) (Fig. 4). These calibrations and our robust phylogenetic trees allow us to investigate the events that shape the divergence of some clades of the order Sepiolida (Figs. 4,  S4).

Sepiolida appeared before the Cretaceous-Paleogene extinction event³⁴, during the middle Mesozoic around 94 Mya (95% HPD = 60.61–130.72). This time frame coincides with the rapid diversification of several oegopsida lineages^15,35. Our molecular estimates also indicate that radiation of Sepiolidae and Sepidariidae occurred around the Cretaceous-Paleogene boundaries and is concurrent with the rapid diversification of modern marine percomorph fishes around the globe, after the extinction of Mesozoic fishes^36,37.

Among the species of Sepiolinae collected in the Mediterranean Sea for this study, only Sepiola robusta Naef, 1912, and Sepiola affinis Naef, 1912 are endemic to the Mediterranean Sea³⁸. The distribution of the other species includes the Mediterranean Sea, North Atlantic Ocean, East Atlantic Ocean, and/or up to the Gulf of Cadiz. The confidence intervals for the split between the Mediterranean-Atlantic and Indo-Pacific lineages, and their diversification, overlap during the early Eocene to the beginning of the Oligocene (Figs. 4 and  S4). This time interval coincides with the end of the Tethys Sea, which separated the Indo-Pacific from the Mediterranean and Atlantic region through the Indian-Mediterranean Seaway^39,40. This separation also influenced the divergence of loliginid clades, coinciding with the split between the Eastern Atlantic plus Mediterranean clade (Loligo, Afrololigo, Alloteuthis) and Indo-Pacific clade (Uroteuthis and Loliolus) (~55 Mya based on Fig. 2 in Anderson and Marian⁴¹).

Our chronogram indicates that the ancestor of Sepiolinae arose prior to the early Eocene around 46 Mya (95% HPD = 25.16–69.49) (Fig. 4), already possessing a bilobed light organ hosting luminescent bacteria (Fig. 3). We estimate that the split between S. affinis and S. intermedia occurred around 2.62 Mya (95% HPD = 0.3–7.4) (Fig. 4) during the end of the Zanclean period, when the Atlantic Ocean refilled the Mediterranean after the Messinian salinity crisis^42,43. While S. affinis is a coastal species with a narrow depth limit, S. intermedia inhabits a wider range of deeper waters. It is possible that two populations of their ancestor, each adapted to a different ecological niche and diverged sympatrically in Mediterranean waters, and, after the speciation, S. intermedia extended its distribution outside the Mediterranean to the Gulf of Cadiz⁴⁴.

We also estimate that the split between H. dispar Rüppell, 1844 and H. hawaiiensis (Berry, 1909) occurred around 2.4 Mya (95% HPD = 0.46–5.88), coinciding with the closure of the Isthmus of Panama around 2.8 Mya⁴⁵. Surveys of these species found H. hawaiiensis in the North Pacific and H. dispar in the North Atlantic Ocean and Mediterranean Sea⁴⁶. A recent speciation event might be the reason for the lack of morphological differences between the two species⁴⁶. Thus, these species may be rendered as cryptic species, a phenomenon increasingly reported in oceanic cephalopods⁴⁷. The sister species of this cryptic species complex, H. dagamensis Robson, 1924, appeared before, around 6 Mya, and is reported with broad distribution in the South Atlantic Ocean off South Africa, the Gulf of Mexico, North Atlantic Ocean between Ireland and Newfoundland in Canada, and the South Pacific Ocean off New Zealand^48,49,50.

The origin of the Heteroteuthis ancestor of H. dispar, H. hawaiiensis, and H. dagamensis can be placed in the Pacific Ocean. After the formation of the Isthmus of Panama, the northern population of Heteroteuthis might have split into H. hawaiiensis in North Pacific and H. dispar in the Atlantic Ocean (from where it also migrated to the Mediterranean Sea). Meanwhile, the formation of the equatorial currents isolated the southern population of Heteroteuthis and gave rise to H. dagamensis. Then, H. dagamensis extended its distribution from the Southern Pacific to the South Atlantic Ocean, the North Atlantic waters, and the Gulf of Mexico. Analysis of molecular species delimitation, however, suggests that H. dagamensis includes cryptic lineages among Atlantic and New Zealand populations³⁰.

While the origin of Heteroteuthis might also be in the Atlantic Ocean, the higher diversity of heteroteuthins in the Pacific (H. hawaiiensis, H. dagamensis, H. ryukyuensis Kubodera, Okutani and Kosuge, 2009, H. nordopacifica Kubodera and Okutani, 2011, and an unknown H. sp. KER (only known from molecular studies⁴⁹)) than at the Atlantic (H. dispar and H. dagamensis), make its origin at the Atlantic less plausible. Moreover, the Atlantic Heteroteuthis were found nested within Heteroteuthinae species from the Pacific, supporting Pacific Ocean origin (Figs. 1, 4).

By sequencing the genomic DNA of sepiolids at low coverage, we recovered complete mitochondrial genomes and nuclear ribosomal genes for most of our collections. Furthermore, mapping reads to the reference genome of E. scolopes allowed us to retrieve additional nuclear-ultraconserved regions. We demonstrate that these nuclear and mitochondrial loci are useful to reconstruct robust phylogenetic trees, especially when the transcriptomes of specimens are difficult to collect, as for sepiolids inhabiting oceanic environments. Finally, our study integrated genomic DNA sequencing with confident morphological identification, which helped to reconstruct the ancestral character of the light organ and its luminescence in sepiolids, and clarify how major lineages have evolved, establishing the existence of distinct Indo-Pacific and Mediterranean-Atlantic subfamilies of Sepiolinae. Our collections and genomically anchored phylogenies will provide a reliable foundation classification of sepiolids for future studies.

Source: Ecology - nature.com