In the molecular evolution of organisms, gene duplication plays a pivotal role in preparing raw materials for evolution [1, 2]. Whereas mutations in single-copy genes are under strong negative selection pressure, redundant gene copies enable a wider spectrum of mutations mostly deleterious. Moreover, some mutations may contribute to sub- or neo-functionalization of copies, making creative evolutionary changes. Gene duplication can occur at the single gene, segmental or chromosomal level, and even whole genome duplications are possible, which may cause the largest scale of divergence of gene functions. In particular, two rounds of whole genome duplication events at initial stages of vertebrate evolution have been important topics for several decades [2–5]. In fish, the ancestor of teleosts (infraclass Teleostei) underwent an additional round of genome duplication (TSD; teleost-specific genome duplication) around 300 million years ago [6–8], which has also attracted the attention of evolutionary biologists. About 48% of extant vertebrate species are teleosts [9, 10], and these live in highly diverse aquatic environments; from tropical coral reefs or rain forests to limit-cooled water under ice in polar regions or highland glacial lakes, from shallow estuaries to abyssal ocean trenches or far offshore open ocean surface, from highly alkalic to acidic waters, and from open waters under the sun to deep into ever dark caves. Thus, the correlation of whole genome duplication and eco-physiological diversification in the teleost lineage has been a topic that has attracted both interest and controversy in evolutionary biology [11–13].

Sensory organs are indispensable for most organisms, of which vision is specifically important for agile animals such as fish. Because of the differential penetration of lights of different wave lengths, diverse light sensitivity may be observed among fish species inhabiting different environments [14]. Vision also physiologically and ontogenetically changes: for example, fish visual sensitivity changes according to habitats shift from freshwater to marine environments [15, 16] and to food shift from carotenoid rich crustacean plankton to blue-green colored pelagic fish [17, 18]. Diversification in light sensitivity at various wave lengths is achieved by evolutionary tuning and differential expression of visual pigment genes. Visual pigment genes have been subjected to studies of gene duplication. In general, fish species have five types of visual pigment genes, which enable effective perception of colors (ultraviolet, blue, green, and red) and dim-light. These subtypes are considered to have arisen in two rounds of whole genome duplications in the ancestor of vertebrates [19]. Regarding TSD, however, the signature remains obscure; one of the two duplicate copies of visual pigment genes that emerged as the result of TSD has been lost in most teleosts. Rather, within each of the five types, the visual pigment gene is often duplicated at the single gene level independently of the teleost lineages [20, 21]. Rhodopsin is the visual pigment protein working for dim light vision [22], and the gene structure has been determined in many animals. The encoding gene (rho) is universally intronless in teleosts, which is believed to be due to a reverse-transcriptional insertion of the spliced messenger RNA transcribed at a distant locus early in the ray-finned fish lineage [23]. The original gene with introns before the insertion is utilized as retinal rhodopsin in tetrapods, but as exorhodopsin in fish [24]. It has been reported that rho is retained in the single copy state in most teleosts [25], but a few exceptional species with duplicate rho are scattered among teleostean tree of life, such as zebrafish, pearl eye, conger, and eel [26–29].

Freshwater eels of the genus Anguilla are among the few catadromous migrating teleost fish, which grow in freshwater areas and descend to the sea for reproduction. Majority of teleostean species live in either marine (55.5%) or freshwater (42.7%) environments, and the other small fraction (1.8%) makes a round trip between the sea and freshwaters for growth and reproduction (diadromy) [30]. Catadromous fish comprise a further small fraction of diadromous fish. The majority of diadromous fish, such as salmons and freshwater gobies, spawn in freshwaters and grow in either marine or freshwater habitats (anadromous or amphidromous). The spawning area of Japanese eel Anguilla japonica Temminck & Schlegel, 1846 has been determined to be located at the southern part of the West Mariana Ridge [31, 32]. After years of living in the freshwater or near-shore habitats in the northeast Asia, Japanese eels metamorphose into silver stage at onset of early maturation and start oceanic migration to reach the above-mentioned open ocean spawning area. Biologging and tracking experiments have revealed that the freshwater eels perform diel vertical migration between upper and lower mesopelagic zones (c.a. 200–1000 m) during oceanic migration (see [33], and references therein). The life cycle of the freshwater eels is thus distinct from other teleosts, which might be associated with their phylogenetic position: the superorder Elopomorpha, as well as Osteoglossomorpha, is considered to have diverged from basal teleosts (Fig. 1) [34, 35].

Ambient light condition of the migrating eels should then change according to their life stages. Carlisle and Denton [16] reported that the wavelength of maximum absorption in visual pigments was different between eels caught in the river and the sea, and further studies have shown that this trait is due to the expression switching of two rho genes in retinal cells [29, 36–38]. Origin of these rho, namely fw-rho for freshwater and ds-rho for deep-sea in this study, is clearly the result of a gene duplication event [29], but when the event occurred remains unsolved. From the deep-branching phylogeny of rhodopsin and related genes in fish, the duplicate genes might be derivatives of TSD, but previous studies have not yielded a clear conclusion [19, 39, 40]. The timing of TSD is estimated to predate the occurrences of Elopomorpha and Osteoglossomorpha, and postdate the divergence between Teleostei and Holostei (gar and bowfin) (Fig. 1). Actually, spotted gar, Lepisosteus oculatus, possesses a single copy of intronless rho in a genomic region on linkage group 5 (LG5), which suggests that this species may have retained the ancestral structure of rho region before TSD [19]. With reference to Osteoglossomorpha, the genome of Asian arowana, Scleropages formosus, was recently sequenced [41], but the structure of its rho region has not been fully investigated.

In this study, we compared eel’s rhodopsin paralogs with those of other teleosts and vertebrates, such as Asian arowana and spotted gar, at the genomic level including synteny organization. The relationship between the evolutionary pathway of rhodopsin genes and migratory ecology is discussed.

We sequenced the genome of Japanese eel, mapped approximately half of the scaffolds assembled to linkage groups. Although the total size of scaffolds was highly consistent with the experimentally estimated genome size [47], it was ~15% larger than the estimate of genome size based on k-mer frequency. Since the genomes sizes of eel species have been recently reconsidered as with their repetitive regions and heterozygous sites [71], further improvements should also be done for the Japanese eel genome. In the present study, we focused on rhodopsin genes in Japanese eel and reconstructed the two loci of rho from the scaffold sequences. The ds-rho and fw-rho were located on different scaffolds, the longer one covering ds-rho was nearly 6 Mb in size, and the shorter one (~400 kb) had at least three paralogs (magi1, prickle2 and atxn7) next to fw-rho. We can thus conclude that the duplication of rho in Japanese eel was neither at the single gene level nor tandem, as in the case of other visual pigment genes, but on a larger scale, and that it was followed by losses of neighboring paralogs. In addition, the scaffold of ds-rho widely corresponded to the region of spotted gar’s LG5, indicating that the genomic region surrounding ds-rho has kept the ancestral structure. Most importantly, we found that the Asian arowana genome had two scaffolds similar to those of ds-rho and fw-rho in Japanese eel, although a copy of rho was lost. Since the both scaffolds widely correspond to spotted gar’s LG5 region, the Asian arowana rho regions are products of a large scale duplication. Eel (Elopomorpha) and arowana (Osteoglossomorpha) are both considered to have been diverged at the early stage of teleost evolution [72]. This suggests a likely scenario for the origin of two regions of rho in Japanese eel and Asian arowana: a large scale duplication occurred in the common ancestor of eel and arowana, and after the divergence of these lineages, both rho genes have been maintained in Japanese eel, whereas a copy of rho has been lost in Asian arowana. An alternative scenario is that large scale duplications have occurred separately in each lineage after the divergences of eel and arowana, respectively, and a copy of rho has been lost in Asian arowana. Based on parsimony, this second scenario is less likely, because it assumes two parallel duplication events. Moreover, it should be noted that TSD is estimated to have occurred shortly before the divergences of eel and arowana [34, 73]. Thus, the second scenario needs further assumptions of loss of either of TSD-derived rho regions and recurrent parallel duplications. In the first scenario on the other hand, the observation can be simply explained by a single large-scale duplication before the divergences of eel and arowana and a loss of rho in the arowana lineage. Considering that the large scaffolds of Asian arowana were mapped on different chromosomes, it is reasonable to conclude that the large-scale duplication was TSD. In this study, we identified ds-rho of Japanese eel on the LG12 according to the mapped linkage markers. In addition, the comparison with the European eel genome suggested that fw-rho might be located on LG6. In the previous study, the LG12 and LG6 of Japanese eel or their backbone regions were estimated to have been built at TSD because of wide correspondences to spotted gar LG5 and medaka’s chromosomes 5 and 7 [49]. Thus, it seems likely that the rho genes in eel have been retained in the descendant regions of TSD-derived chromosomes.

Scaffold 435 of fw-rho in Japanese eel is short and its correspondence to Asian arowana’s scaffold 11 remains obscure, raising the additional possibility that fw-rho in Japanese eel has been generated by another duplication (e.g., lineage-specific segmental duplication). Although this third scenario also assumes the loss and a recurrent duplication of the rho region after TSD as in the second scenario, multiple losses and duplications may be possible on a smaller scale, such as the tandem duplications reported in other visual pigment genes. Therefore, we further examined the microsynteny around rho locus among Japanese eel and other teleosts. In spotted gar as an outgroup, rho was located within the region from lrig1 to thoc7 in LG5, and the synteny was also observed in the corresponding regions of Japanese eel and Asian arowana, strongly suggesting that the gene order is the ancestral form in teleosts. Considering that paralogs of magi1, prickle2 and atxn7 are present in the fw-rho region of Japanese eel, at least eight genes (magi1, ds-rho, adamst9, prickle2, pphln1, slc2a9l1, psmd6, and atxn7) as a template should have been doubled in the duplication event, followed by losses of several paralogs. Apparently, this is also the case in Asian arowana. Here, it is important to note that the pattern of gene losses was the same between Japanese eel and Asian arowana, except for that a copy of rho is lost in Asian arowana. Many studies have shown that either of the paralogs generated by duplication will be rapidly lost, while a number of paralogs might diverge functionally to each other [74–76]. The pattern of gene losses is a good signature for estimating the evolutionary scenario of duplicate genomic regions [77], because the genes lost once in the template region would never be restored by duplication again. In this study, for example, we can say that the fw-rho region with pphln1 in Japanese eel should have been generated before the pphln1 paralog was lost in the ds-rho region, and so too with those in Asian arowana. In particular, four genes, adamst9, pphln1, slc2a9l1, and psmd6, were lost in the same pattern between Japanese eel and Asian arowana. Assuming that the duplications occurred separately in each of eel and arowana lineages, these coincidences are unlikely. Moreover, in all the teleost genomes examined, we did not find any cases inconsistent with this pattern. Thus, these results suggest that the gene orders in two rho regions observed in Japanese eel had already emerged in the ancestor of eel, arowana, and clupeocephalan (Otocephala + Euteleostei) lineages (Fig. 5). Whereas zebrafish and Mexican tetra have two copies of rho, all the euteleosts examined have a single copy of rho, but the gene loss was in the opposite side from that in Asian arowana. This is also explained by our hypothesis that two rho genes were still present at the divergence of teleostean lineages. Thus, we propose that the losses of rho on opposite regions occurred independently in arowana and clupeocephalan lineages after their divergence; otherwise it would be difficult to explain the loss pattern of rho and adamst9 paralogs between magi1 and prickle2. Regarding two rho copies of zebrafish, phylogenetic analysis has recently suggested that the duplication occurred early in a teleost lineage, and either of the paralogs has evolved under positive selection [39]. According to our hypothesis, considering that the zebrafish lineage had branched at the base of Clupeocephala, it is possible that two rho copies predating the divergence of teleostean lineages have remained till the present in zebrafish and lost in the euteleostean ancestor. This possibility may be further tested by comparison to two rho copies in Mexican tetra in the same superorder (Otophysi) as zebrafish.

Two synteny clusters of the rho region among 11 teleosts in trees based on concatenated sequences (Figs. 4c, d) were in accord with the first scenario, because mutual monophyly of these synteny clusters can be traced back to the teleostean root. Recovery of outline of the teleostean interrelationships (clupeocephalan and euteleostean clades) in both synteny clusters with high branch supports indicates that the trees obtained from concatenated sequences are close to the true tree. The results from concatenated sequences, however, do not necessarily support the first scenario, because the two paralogous rho regions are merely vehicles of the gene of interest. Clustering of the eel rho copies relative to other teleosts in trees solely based on rho sequences regardless of exclusion of sites with higher non-synonymous substitutions (Figs. 4a, b and Additional File 10) is possible upon delayed resolution of tetrasomy [78] after TSD. Tetrasomy of chromosomal segments homogenizes homeolog sequences through recombination. Phylogenetic affinity of eel’s two rho copies could arise, if the rho region was tetrasomic and the resolution (evolution of paralogy) was delayed particularly in the eel lineage. Atxn7-based trees also show clusters of eel and/or arowana paralogs (Additional file 9). Allotetraploidy (no tetrasomy predicted) at TSD is suggested [79], but tetrasomy can occur even after allotetraploidization [80]. Another possibility for this skewed tree shape arises upon the second scenario with recurrent gene losses/duplications.

Nevertheless, dense taxon sampling of rho sequences [39], magi1 and prickle2 sequences in our study (Additional file 9) yielded similar results with those based on concatenated sequences. Differences of phylogenetic resolving powers in different sizes of sequence datasets (either or both dimensions of length and number of OTUs) would be a source of the incongruences. Our results suggest affinity of ds-rho and fw-rho with synteny 1 and 2 clusters respectively, and thus the origin of ds- and fw-rho genes can be traced back to the basal teleostean divergence near TSD. Thus, phylogenetic analysis supports the hypothesis that rho paralogs in Japanese eel have originated from TSD.

Finally, our results might shed light on the relationship between TSD and evolution of elopomorph species. Regarding sexual maturation, hormone receptor paralogs built at TSD, which are present in most other teleosts, have been examined also in eels [77, 81–83]. Moreover, genome-wide studies have shown that eels have the complete set of Hox gene clusters doubled at TSD [60, 84], while other teleosts lack some of the genes. Although such duplicate genes are likely involved in eels’ life cycle, the genes are largely conserved among teleosts and might not be directly correlated with eels’ unique adaptation to the freshwater and marine environments. Since the two rhodopsins in eel have been tuned for these environments, respectively [29, 36–38], our finding provides the first evidence about the functional differentiation of TSD-derived paralogs with special reference to the eel-specific ecological traits. The species which had deviated earlier in elopomorph divergence are distributed from estuary (e.g., tarpon) to deep sea (e.g., gissu). After their divergence, the lineages of eels and relatives (conger, gulper eel, etc.) are inferred to have diverged [73, 85, 86]. These studies imply that the aquatic habitats of elopomorph species have been often changed. Therefore, the present study raises the interest in the correlation between such a habitat alteration and the functions of rhodopsin in elopomorph species. For example, whether two copies of rho are maintained or either is lost in these species would be a simple but important topic. Japanese Conger conger myriaster has two copies of rho orthologous to those of Japanese eel, respectively [28], hence these genes are also derived from TSD as shown in this study. Previous studies showed functional differentiation of these rho copies in Japanese conger but in a different way from Japanese eel [28, 70, 87], suggesting multiple occurrence of the functional differentiation of rho paralogs by minor mutations in elopomorph fish retaining paralogs in their genome. In Anguilla, such as Japanese eel and European eel, the ancestral forms of rho may have been both deep-sea types, one of which has been mutated toward a freshwater-type, because a recent study suggested that these species evolved from a deep-sea habitat [88]. Thus, maintenance and sub-functionalization of the rho paralogs might be correlated with the wide distribution of elopomorph species from freshwater to the deep sea. It is difficult to answer in this study due to lack of other genome data on how often the functional changes of rho paralogs occurred in Elopomorpha. However, it should be stressed that maintenance of TSD-derived rho paralogs may have finally contributed to the establishment of eel’s life cycle, particularly the habitat alteration by long-distance migration between distinct aquatic environments. Further research based on genome-wide comparison may reveal whether or not this is a rare case with regard to the evolutionary impact of TSD.

Two copies of eel rhodopsin gene were generated most likely at TSD. This is the first finding that TSD (~300 million years ago) has remotely lead to the gene innovation involved in eel-specific migratory life cycle. Further researches may provide insight into correlations between the functional differentiation of TSD-derived paralogs and the diversification of elopomorph species.