Skip to main content
Download PDF

Expression variations in ectodysplasin-A gene (eda) may contribute to morphological divergence of scales in haplochromine cichlids

BMC Ecology and Evolution volume 22, Article number: 28 (2022) Cite this article

Abstract

Background

Elasmoid scales are one of the most common dermal appendages and can be found in almost all species of bony fish differing greatly in their shape. Whilst the genetic underpinnings behind elasmoid scale development have been investigated, not much is known about the mechanisms involved in moulding of scales. To investigate the links between gene expression differences and morphological divergence, we inferred shape variation of scales from two different areas of the body (anterior and posterior) stemming from ten haplochromine cichlid species from different origins (Lake Tanganyika, Lake Malawi, Lake Victoria and riverine). Additionally, we investigated transcriptional differences of a set of genes known to be involved in scale development and morphogenesis in fish.

Results

We found that scales from the anterior and posterior part of the body strongly differ in their overall shape, and a separate look on scales from each body part revealed similar trajectories of shape differences considering the lake origin of single investigated species. Above all, nine as well as 11 out of 16 target genes showed expression differences between the lakes for the anterior and posterior dataset, respectively. Whereas in posterior scales four genes (dlx5, eda, rankl and shh) revealed significant correlations between expression and morphological differentiation, in anterior scales only one gene (eda) showed such a correlation. Furthermore, eda displayed the most significant expression difference between species of Lake Tanganyika and species of the other two younger lakes. Finally, we found genetic differences in downstream regions of eda gene (e.g., in the eda-tnfsf13b inter-genic region) that are associated with observed expression differences. This is reminiscent of a genetic difference in the eda-tnfsf13b inter-genic region which leads to gain or loss of armour plates in stickleback.

Conclusion

These findings provide evidence for cross-species transcriptional differences of an important morphogenetic factor, eda, which is involved in formation of ectodermal appendages. These expression differences appeared to be associated with morphological differences observed in the scales of haplochromine cichlids indicating potential role of eda mediated signal in divergent scale morphogenesis in fish.

Peer Review reports

Background

Cichlids pose a great a model system for evolutionary biology, as they include some of the most striking examples of explosive adaptive radiation and many aspects of their life history as well as their behaviour, coloration and feeding morphologies are well studied [1,2,3]. One of the most remarkable features is their repeated evolution of parallel eco-morphologies, especially across the radiations of the three East African Great Lakes, Lake Tanganyika (LT), Lake Malawi (LM) and Lake Victoria (LV) [4, 5]. These ecological adaptations are also the focus of many studies, as they promise the opportunity to shed light on different molecular mechanisms underlying repeated evolution and diversification [6, 7]. Regarding skeletal morphogenesis in particular the evolution of their jaws and their phenotypic plasticity are topics of ongoing research [7,8,9,10,11,12]. However, while the adaptive value of some of the investigated structures (e.g., feeding apparatus) can be more easily connected to certain ecological specializations [5, 13], this is not so obvious in others, such as scales.

Fish scales come in a vast array of different shapes and forms. As a part of the dermal skeleton, which amongst other structures also includes teeth, odontodes, spines and fin rays, these postcranial derivates evolved into morphologically and histologically diverse structures in Actinopterygii [14, 15]. Elasmoid scales, found in most of teleost species, form in the dermal mesenchyme and are mainly used for protection and hypothetically for hydrodynamic modifications [14, 16, 17]. While the elasmoid scales form relatively late in ontogeny and can take diverse forms, they share a composition consisting of three tissues, with elasmodin as the basal component formed in a characteristic plywood-like structure [15, 16]. Scale development, mostly studied in zebrafish, has been found to be orchestrated by several well-known pathways, including Hh, Fgf and Eda [16, 18,19,20], which are known to be also involved in the appendage formation across several vertebrate groups [21]. Mutations and allele variations in the Eda/Edar pathway, for example, have been linked to fish fin, scale and armour plate development as well as human and mouse hair and teeth growth [19, 22, 23]. Nevertheless, besides a recent extensive comparison of the scale morphology across Lake Tanganyika cichlids [24], as well as a genetic study of scale shapes in two closely related Lake Malawi cichlids, which tied FgF signalling to scale shape variation [20], not much is known about the molecular mechanisms shaping the elasmoid scale.

In this study, we investigate the morphological differences in the anterior and posterior scales of ten haplochromine cichlid fish species from three Great East African Lakes, i.e., Lake Tanganyika (LT), Lake Malawi (LM) and Lake Victoria (LV) as well as a riverine haplochromine cichlid species. After identification of a stably expressed reference gene, we also investigate transcriptional differences of a set of genes known to be involved in scale development and morphogenesis in fish. Finally, we tried to find links between the gene expression differences and morphological divergence in both anterior and posterior scales. Our results provide cross-species expression comparisons of scale related genes in haplochromine cichlids and implicate expression differences by which formation of distinct scale morphologies might be determined.

Methods

Fish husbandry and sampling

Ten haplochromine cichlid species; three species from Lake Tanganyika, four species from Lake Malawi, two species from Lake Victoria, and one riverine haplochromine species, were selected for this study (Fig. 1A). For each species, the first generation siblings of wild-caught fish from the aquarium trade were bred and raised in standardized tanks and rearing conditions with the same diet (Spirulina flakes) until they displayed mating behaviour. Between five to 11 adult males per species were sampled for morphological analysis and four adult males were sampled for gene expression investigation. The sampled fish species were sacrificed by euthanization with 0.5 g MS-222/litre of water, and five anterior and posterior scales from left side of the body were removed for morphological analysis (Fig. 1B), whereas similar numbers of scales (together with their attached covering epidermis) were taken from both sides and all anterior or posterior scales from each fish were pooled for gene expression analysis.

Fig. 1
figure 1

The haplochromine cichlid species and descriptions of the scale samples. a A simplified phylogenetic relatedness of the East African haplochromine cichlid species used in this study. b Positions of the anterior and posterior scales used in this study and landmarks used for the geometric morphometric analyses in both anterior (left) and posterior (right) scales shown as an example for Petrochromis famula. Blue dots represent major landmarks, white dots semi-landmarks and 1 mm scale bars are given below the images. c Principal component analysis (PCA) plots clearly separate scales from the anterior and posterior part of the body. Additional, warped outline drawings illustrate major shape changes along the axis (red) compared to the overall mean shape (grey)

Full size image

Morphological analysis

To infer shape differences of scales from divergent African cichlids from different lakes a 2D geometric morphometric framework was deployed. Due to major morphological differences of the scales they were separately investigated for the anterior and posterior part of the body (Fig. 1 B and C). Standardized images of scales were taken with a KEYENCE VHX-5000 digital microscope (KEYENCE Germany GmbH). 84 adult specimens from 10 cichlid species inhabiting the three major rift lakes which were reared under standardized aquarium conditions (Astatotilapia burtoni = 7; Neochromis omnicaeruleus = 10; Petrochromis famula = 11, P. polyodon = 7; Paralabidochromis sauvage = 5; Simochromis diagramma = 11; Sciaenochromis fryeri = 5; Tropheops tropheops = 9; Labeotropheus trewavasae = 9; Mz: Maylandia zebra = 10) were included for the geometric morphometric analyses. For each individual six scale replicates from the anterior and posterior part of the body were probed (Fig. 1B), leading to a total of 1.008 investigated scales. After randomizing pictures in tpsUtil v.1.6 (available at http://life.bio. sunysb.edu/morph/soft-utility.html), landmark digitization was conducted on a set of seven fixed landmarks and 14 semi-landmarks (see Fig. 1b for positions) in tpsDig v.2.26 (available at http://life.bio.sunysb.edu/ morph/soft-utility.html). To ensure consistency, this step was conducted by a single investigator. Generalized Procrustes superimposition [25] was performed in tpsRelw v.1.65 (available at http://life.bio.sunysb.edu/morph/soft- utility.html) and aligned landmark configurations were exported for further analysis in MorphoJ v.1.06 [26]. In MorphoJ, single observations obtained from the six replicates were averaged to get the mean shape for each landmark. A Principal Component analysis (PCA) was applied to infer variation in morphospace among scale position (anterior vs. posterior), single specimen, and species. Subsequent analyses were based on separated datasets for anterior and posterior scale landmark setting, whereas PC-scores were exported for linear discriminant function analyses (LDA) in PAST v.4.1 [27]. To reduce the number of variables and control for putative over-separation of groups [28], only the first four principal components were used for the LDA. PCA and LDA plots were visualized in R v3.1.2 [29].

RNA isolation and cDNA synthesis

As mentioned in the section above, for isolating the total RNA ten anterior and posterior scales from each fish were pooled in a single tube containing 0.25 mL of a tissue lysis buffer from Reliaprep RNA tissue miniprep system (Promega, #Z6111, USA) as well as one 1.4 mm ceramic bead to crush the scales. The scales were homogenized using a FastPrep-24 Instrument (MP Biomedicals, Santa Ana, CA, USA) and total RNA was extracted following the instructions provided by the manufacturer (adjusted protocol for small amounts of fibrous tissue). In summary, the instructions proceed with mixing of the lysis buffer and homogenized scales with isopropanol and centrifuging the entire mix through a column provided by the kit, followed by several RNA washing steps and a final DNase treatment. The RNAs were quantified by a Nanophotometer (IMPLEN GmbH, Munich, Germany) and their quality was checked with RNA ScreenTapes on an Agilent 2200 TapeStation (Agilent Technologies). Next, the RNA samples with a RNA integrity number (RIN) above six were applied to first strand cDNA synthesis using 300 ng of RNA and High Capacity cDNA Reverse Transcription kit (Applied Biosystems). The synthesized cDNA from each RNA sample was diluted 1:5 times in nuclease-free water to conduct qPCR.

Gene selection, designing primers and binding site predictions

We selected eight candidate reference genes which have been frequently used in different studies of Haplochromine cichlids and have shown high expression levels in various connective tissues including skeletal tissues [10, 30,31,32,33,34,35]. Furthermore, we chose 16 target candidate genes, which are implicated in scale development and morphogenesis (Table 1). The primers were designed at conserved sequence of coding regions using already available transcriptome/genome data of East African haplochromine species on www.ncbi.nlm.nih.gov/sra, including Astatotilapia burtoni (SRX4523155), Labeotropheus trewavasae (SRX6432658), Metriaclima zebra (SRX8892877), Neochromis omnicaeruleus (SRX8567938), Paralabidochromis sauvagei (SRX8892920), Petrochromis famula (SRX6445829), Petrochromis polyodon (ERX3501455), Sciaenochromis benthicola (ERX1818621), Simochromis diagramma (SRX8567938), and Tropheops tropheops (ERX659441) as well as two more distant species from different African cichlid tribes (Oreochromis niloticus and Neolamprologus brichardi) [7, 36,37,38,39]. The sequences from all the species were imported to CLC Genomic Workbench, version 7.5 (CLC Bio, Aarhus, Denmark), and after alignment, the exon/exon junctions were specified using the Astatotilapia burtoni annotated genome in the Ensembl database (http://www.ensembl.org) [40]. The primers were designed spanning exon junctions and a short amplicon size (< 250 bp) as recommended to be optimal for qPCR quantification [41]. The primers were designed and assessed through Primer Express 3.0 (Applied Biosystems, CA, USA) and OligoAnalyzer 3.1 (Integrated DNA Technology) to minimize the occurrence of dimerization and secondary structures.

Table 1 Selected target genes involved in the development and/or morphogenesis of scales in teleost fish
Full size table

We retrieved downstream sequences (3’UTR and inter-genic region) of eda gene for all the species in this study from European Nucleotide Archive (ENA) and Sequence Read Archive (SRA) in order to identify changes in potential binding sites. To do this, we used genomic sequences of the haplochromine species; A. burtoni (GCA_000239415.1), P. famula (GCA_015108095.1), N. omnicaeruleus (SRR12700904), P. polyodon (GCA_015103895.1), S. fryeri (ERX1818621), S. diagramma (GCA_900408965.1), M. zebra (GCA_000238955.1), T. tropheops (SAMEA2661272), L. trewavasae (SAMN12216683), and P. sauvagei (GCA_018403495.1). Next we identify the 3’UTR and inter-genic region of eda genes using the annotated genome of A. burtoni from Ensembl and aligned them using CLC Genomic Workbench. The different sequence motifs were identified and screened for potential TF binding sites using STAMP [42] and the PWMs obtained from the TRANSFAC database [43].

qPCR and data analysis

The qPCR reactions were generated using Maxima SYBR Green/ROX qPCR Master Mix (2X) (Thermo Fisher Scientific, Germany) and the amplifications were conducted on ABI 7500 real-time PCR System (Applied Biosystems). The qPCR setups followed the recommended optimal sample maximization method [44]. The qPCR program, dissociation step and calculation of primer efficiencies were performed as described in our previous study [45] (Additional file 1).

Three different algorithms were applied to validate the most stable reference genes; BestKeeper [46], NormFinder [47] and geNorm [48]. The Cq value of the most stable reference gene was used as normalization factor (Cq reference) to calculate ΔCq of each target gene (ΔCq target = Cq target – Cq reference). The lowest expressed sample in each expression comparison was used as a calibrator sample and rest of the samples were subtracted from its ΔCq value to calculate ΔΔCq values (ΔCq target – ΔCq calibrator). Relative expression quantities (RQ) were calculated through E−ΔΔCq [49]. In order to perform statistical analysis, fold differences (FD) were calculated by transformation of RQ values to logarithmic values [50]. The significant expression differences were calculated using ANOVA statistical tests, followed by Tukey’s HSD post hoc tests. The correlations between gene expression and a morphometric parameter (canonical variate 1) were calculated through Pearson correlation coefficients (r) for each gene using R.

Results

Divergence in scale morphology

The principal component analysis (PCA) revealed a clear separation in overall average individual shape between anterior and posterior scales (Fig. 1C). PC1 and PC2 explained 78.3% and 8.6% of the total shape variation, respectively. Generally, on the first axis anterior scales are anterior-posteriorly more compressed compared to the posterior body part (see deformation grids in Fig. 1C). Along the second PC axis changes can be observed in the shape of the posterior scale field (narrow vs. wide), as well as in the lateral edges of the anterior scale field (edges vs. round).

While comparing different species, large variation in overall shape can be observed in the dataset which is restricted to anterior scales only. Petrochromis famula, Maylandia zebra and Sciaenochromis fryeri occupy large parts of the morphospace and overall, less intraspecific variation can be observed in other species (Fig. 2A). In the anterior dataset changes along the PC1 explain 47.1% of total variation, and mainly affect the circularity of the overall shape (i.e., that scales get more compressed towards positive values). Changes along the second PC, which explains 19.6% of the total shape variation, affect the posterior scale field (compression vs. expansion). Compared to anterior scales, less intraspecific shape variation can be observed in posterior scales, whereas PC1 explains 54.8% and PC2 18.9% of the total variation, respectively (Fig. 2B). PC1 separates two major clusters (S. diagramma + P. famula vs. rest) whereas changes along the axis mainly contribute to a dorso-ventrally versus anterior–posterior compression of the scale and the roundness of the anterior scale field. Along PC2 shape changes affect the expansion (or compression) of the anterior and posterior scale fields. Generally, the PCA only poorly resolves the lake (or phylogenetic) origin of the single species for both the anterior and posterior dataset.

Fig. 2
figure 2

Morphospace of investigated scales from different species and body parts. Principal component analysis based on average shape of scales collected for the anterior a and posterior b part of the body and respective shape differences along the axis (grey: overall mean shape; red: shape change). Linear discriminant function analysis based on the first four PC-scores for anterior c and posterior d scales. All data points represent mean shapes obtained from 6 individually collected scales and shapes represent different lake origins. A.b. Astatotilapia burtoni, N.o. Neochromis omnicaeruleus, P.f. Petrochromis famula, P.p. P. polyodon, P.s. Paralabidochromis sauvage, S.d. Simochromis diagramma, S.f. Sciaenochromis fryeri, T.t. Tropheops tropheops, L.t. Labeotropheus trewavasae, M.z. Maylandia zebra

Full size image

The linear discriminant function analysis (LDA) of anterior as well as posterior dataset correctly classified 77.38% (jackknifed: 67.86%) and 72.62% (jackknifed: 65.48%) of species (Fig. 2C, D). The first axis explains 83.66% and 76.00% variance of the overall shape variability for the anterior and posterior dataset, respectively. In the anterior dataset, the first LD-axis separates three major clusters made up of samples from Lake Tanganyika, the riverine Astatotilapia burtoni and a joint Victoria-Malawi cluster. Similar results were obtained for the posterior dataset, whereas along the first axis the separation between the riverine A. burtoni and the Victoria-Malawi cluster is less prominent. Along the second axis, which explains 10.04% and 13.2% of the variance of the overall shape variability for the anterior and posterior dataset, respectively, mainly interspecific and intraspecific variation is portrayed. Overall, for the anterior dataset, 83.33% (jackknifed: 75%) of species were correctly classified according to the lake origin, whereas single classification scores reached values of 100% (jackknifed: 85.71%) for A. burtoni, as well as 66.67% (jackknifed: 57.58%) for Malawi, 86.67% (jackknifed: 86.67%) for Victoria and 96.55% (jackknifed: 86.21%). In total, for the posterior dataset, 77.38% (jackknifed: 70.24%) of the individuals were correctly assigned to the lake origin, with 85.71% (jackknifed: 71.43%) for A. burtoni, as well as 66.67% (jackknifed: 54.55%) for Malawi, 53.33% (jackknifed: 46.67%) for Victoria and 100% (jackknifed: 100%) for individuals from Tanganyika.

Validation of stable reference genes

To quantify the expression levels of the selected target genes, the validation of stable reference gene(s) with least variation in expression across the anterior and posterior scales of different species is a necessary step [51]. The eight candidates were selected from frequently used reference genes in studies of different tissues in East African cichlids [10, 30,31,32,33,34,35]. The candidate reference genes showed variable expression levels in the scales, and from highest to lowest expressed were respectively; actb1, hsp90a, rps11, rps18, hprt1,gapdh, elf1a and tbp. Interestingly, in both anterior and posterior scales, all the three software ranked actb1 as the most stable reference gene with lowest expression variation across the cichlid species in this study (Table 2). Thus, we used the Cq value of actb1 as normalization factor (NF) in each sample for quantification of relative expression analyses of the target genes.

Table 2 Ranking of reference genes in anterior and posterior scales across all of the haplochromine species used in this study
Full size table

Gene expression differences between anterior and posterior scales

The relative expression levels of 16 candidate target genes, bmp4, col1a2, ctsk, dlx5, eda, edar, fgf20, fgfr1, mmp2, mmp9, opg, rankl, runx2a, sema4d, shh and sp7, were compared between the anterior and posterior scales in each of the haplochromine species (Fig. 3). Some of these genes such as bmp4, col1a2, rankl and sp7, showed almost no expression difference between the anterior and posterior scales. Moreover, none of the target genes showed consistent expression difference across all the species. These indicate potential involvement of various genes in morphological divergence between the anterior and posterior scales. However, two genes, ctsk and shh exhibited expression difference between the anterior and posterior scales in most of the species (Fig. 3). The directions of expression differences between the anterior and posterior scales for ctsk and shh were variable depending on the species. Interestingly all the three species form LT showed higher expression in the anterior scale for shh, whereas the all the species from LM and LV showed tendency for opposite pattern with increased posterior scale expression. These findings suggest potential role of ctsk and shh in morphological divergence of the scales along anterior–posterior axis.

Fig. 3
figure 3

The anterior versus posterior scales expression differences of the candidate target genes in haplochromine cichlids from three East African lakes. Comparisons of relative expression levels between anterior versus posterior scales for 16 candidate target genes in different lakes in East Africa at young adult stage. Significant differences between the are indicated by red asterisks (*P < 0.05; **P < 0.01). See Fig. 1A for corresponding species abbreviations

Full size image

Gene expression differences between lakes in anterior and posterior scales

Next, we compared the expression levels of the target genes between the lakes in the anterior or posterior scales by considering all the species from each lake as one group (Fig. 4). In the anterior scales, nine out of 16 target genes showed expression differences between the lakes including bmp4, ctsk, eda, edar, mmp2, opg, rankl, shh and sp7. Most of these differences were between LT and the other lakes, and 4 genes, ctsk, mmp2, opg and rankl showed higher expression in LT species, while 3 genes, bmp4, eda and shh showed lower expression LT species. Furthermore, 2 genes, ctsk and eda displayed the strongest expression differences between the lakes in opposite patterns suggesting their role in morphological divergence of the anterior scales across the lakes (Fig. 4). In the posterior scales, 11 out of 16 target genes showed expression differences between the lakes including bmp4, col1a2, ctsk, eda, edar, fgf20, fgfr1, mmp2, opg, rankl, shh and sp7. Again, most of these differences in the posterior scales were between LT and the other lakes, and five genes, col1a2, ctsk, mmp2, opg and rankl showed higher expression in LT species, while three genes, eda, fgf20 and shh showed lower expression in LT species. In addition, four genes, fgf20, rankl, eda and shh displayed the strongest expression differences between the lakes in opposite patterns (eda and rankl higher in LT, and fgf20 and shh lower in LT) suggesting their role in morphological divergence of the posterior scales across the lakes (Fig. 4). In general, more genes with stronger expression differences between the lakes were observed the posterior scales. Several genes such as ctsk, eda, edar, mmp2, opg, rankl and shh appeared to have similar patterns of expression differences between the lakes in both anterior and posterior scales. Importantly, we found only one gene, eda, to have strong differences between the lakes in both anterior and posterior scales indicating its potentially crucial role in morphological divergence of the scales across the lakes.

Fig. 4
figure 4

The lake-based expression differences of the candidate target genes in haplochromine cichlids in this study. Comparisons of relative expression levels between the lakes, when all species of each lake were combined, within anterior or posterior scales for 16 candidate target genes. Significant differences between the are indicated by red asterisks (*P < 0.05; **P < 0.01; ***P < 0.001). See Fig. 1A for corresponding species abbreviations

Full size image

Correlation analyses between gene expression and morphological divergence in scales

We analysed the correlation between expression of the genes and canonical variate 1 in the anterior or posterior scales across the species. Only one gene, eda, showed significant correlation in the anterior scales (Fig. 5), whereas, in the posterior scales four genes including dlx5, eda, rankl and shh displayed significant correlations between expression and morphological differences (Fig. 6). Among these genes eda exhibited the strongest correlation in the posterior scales. However, the correlation patterns differed between the genes in the posterior scales, i.e., eda and shh showed positive while dlx5 and rankl had negative correlations with the morphological changes based canonical variate 1. Therefore, again only one gene, eda, showed significant correlation between its expression and the morphological differences in both scales indicating its potential role in divergent scale morphogenesis in the cichlid species. The opposite correlation patterns in the posterior scales might also indicate inhibitory regulatory connections between the genes.

Fig. 5
figure 5

Correlation analyses of candidate target gene expressions and the anterior scale morphololgical divergence across the haplochromine species. A Pearson correlation coefficient (r) was used to assess the similarity between differences in expression level of the target genes and the major canical variate in the anterior scales across all species. See Fig. 1A for corresponding species abbreviations

Full size image
Fig. 6
figure 6

Correlation analyses of candidate target gene expressions and the posterior scale morphololgical divergence across the haplochromine species. A Pearson correlation coefficient (r) was used to assess the similarity between differences in expression level of the target genes and the major canical variate in the posterior scales across all species. See Fig. 1A for corresponding species abbreviations

Full size image

Genetic differences in non-coding sequences of eda gene

Finally, we were interested to investigate genetic differences in available regulatory sequences of eda gene including 5´UTR, 3´UTR and short but conserved inter-genic region between eda and tnfsf13b (immediate downstream gene) across the species. Interestingly, we found two genetic differences (mutations/deletions) in 3´UTR and one in eda-tnfsf13b inter-genic region to differ between LT species versus LM and LV species (Table 3). Next, we parsed the short sequence regions containing the mutations/deletions against transcription factor binding site (TFBS) databases. We found that the two changes in 3´UTR seem to lead to gaining TFBS for transcription factors Mef2 and Tcf1 in the LM and LV species, whereas the changes in the inter-genic region led to gaining TFBS for Lef1 transcription factor in the LT species (Table 3). Importantly, the riverine species A. burtoni appeared to have intermediate genetic changes meaning that for the two changes in 3’UTR it showed a deletion similar to the LT species but a mutation similar to the LV and LM species. Also, for the inter-genic change, A. burtoni showed an intermediate mutation between the LT and the other species from LM and LV, however, this mutation showed no gain of TFBS (similar to the LM and LV species). Taken together, these genetic changes showed similarity with differences in gene expression and scale morphology, where the LT species clustered different from LM and LV species and the riverine species (A. burtoni) showed intermediate differences. This suggests that the identified genetic changes might be the underlying factors for divergent eda expression as well as differences in the scale morphology.

Table 3 Identified genetic differences in regulatory sequences of eda gene and predicted binding sites for potential upstream regulators
Full size table

Discussion

As river-adapted haplochromine cichlids repeatedly seeded adaptive radiations in several East African lakes, cichlid fishes recurrently adapted to corresponding trophic niches. Thereby, the adaptive value of traits is often mirrored by morphological shape parallelism and concomitant similar lifestyles which result from parallel evolution [4, 5]. Hence, fishes from the cichlid species flocks in various African lakes comprise an exciting model to conduct comparative morphological and molecular studies. While most previous studies focused on bony elements that can easily be linked particular trophic niches and divergent natural selection as driver of diversification in cichlids (e.g., [52]), other skeletal structures such as scales might show less obvious adaptive trajectories.

Above all, between the three East African Great Lakes, the haplochromine cichlids are especially interesting, as they share common ancestry and comprise the Tropheini at LT and the entire the LV and LM haplochromine radiations [53, 54]. As the lakes have all very different geological histories, with Lake Tanganyika being the oldest [55], LM the intermediate [56] and LV the youngest of the three [57], they also depict three extensive radiations at different time points. Thus, depending on the evolutionary age of the different lakes, species (and their morphologies) had more or less time to diverge, despite sharing parallels. The more time passes, much more elaborated predator–prey and host-parasite relationships can evolve. This is manifested by unique ecological and behavioural features, particularly in the oldest of the three lakes, LT, which contains coocoo-catfish species showing brood parasitism [58], dwarfed gastropod shell breeders [59], putative cleaning behaviour [60], or highly elaborated scale eaters [61] (but also see Genyochromis mento from Lake Malawi). Particularly the latter case, scale eating, could have influenced the co-evolution of scale morphology in host species. Lake Tanganyika’s scale eaters (i.e., Perissodus; Perissodini) show different degrees of specialization, whereas only the shallow water species, Perrisodus microlepis and P. straeleni, feed almost exclusively on fish scales while other species are not that specialized [4, 61]. The most common prey species of P. microlepis are members of the Tropheini and Eretmodini [62]. P. straeleni seems to be less specialised to certain prey items, but Tropheini scales still make up a major part of the gut contents [63]. Based on our dataset (Tropheini only) it remains speculative to assume that scale-related gene expression and the concomitant morphology might reflect an adaptation to reduce the risk of scale predation. Future studies, including more early branching (non-modern) haplochromine cichlids (e.g., Pseudocrenilabrus, Thoracochromis, Astatoreochromis) or other Tropheini with different lifestyles (e.g., Ctenochromis) will be necessary to establish stronger links between scale morphology and this unique predation pressure.

Nonetheless, understanding which genetic mechanisms underlie the scale morphology might be the key to understand how such similar and/or divergent eco-morphologies evolved. Perhaps the most striking finding of our study was the highly significant differential expression of eda between LT species versus the species from the younger lakes (LM and LV species) in both anterior and posterior scales (Fig. 4). Interestingly, the eda expression in the riverine species (A. burtoni), which recently has been shown to be closer related to the LM and LV species flocks as a part of the sister group to the LT Tropheini [64], was at intermediate level between LT and the species from LM and LV in both anterior and posterior scales, reflecting the phylogeny. Moreover, the expression patterns of eda in both scales were highly correlated with morphological divergence across the species in this study. Ectodysplasin A (eda) encodes a member of the tumor necrosis factor family and mediates a signal conserved across vertebrates which is essential for morphogenesis of ectodermal appendages, such as scale, hair and feathers [22]. The eda signal is mediated through its receptor (encoded by edar) and initiated upon binding of eda to edar on the surface of a target cell [22]. In human, mutations in components of eda signal can cause hypohidrotic ectodermal dysplasia (HED) which is characterized by reduction and abnormal teeth morphology, absence or reduction of sweating glands and hair [65]. Similarly, impaired eda signal in zebrafish and medaka can lead to reduction in the number of scales and teeth [18, 19]. In sculpin (Cottus) fishes, genetic changes in the receptor gene (edar) has been found to be associated with morphological variations in body prickles (calcified spicules embedded in the skin), which are homologous structures to fish scales [66]. A later study in a highly derived order of teleosts, Tetraodontiformes, which includes ocean sunfishes, triggerfishes and pufferfishes, also showed the importance of eda signaling pathway in developmental formation and morphological variations of dermal spines (an extreme scale derivative) [67]. In stickleback, a mutation within an inter-genic region between eda and tnfsf13b genes leads to changes in transcriptional responsiveness of eda to its upstream Wnt signaling pathway and consequently impairment of armor plate formation [68].

In this study, we also found genetic changes in 3’UTR and the inter-genic region between eda and tnfsf13b genes that could explain the differences in eda expression across the Haplochromine cichlids (Table 3). The genetic changes resulted gain or loss of motifs which were predicted to be binding sites for transcription factors encoded by mef2, tcf1 and lef1 genes. These changes always discriminated the LT species from the species from LM and LV, whereas the riverine species had changes which could be considered an intermediate to both groups. Interestingly, all of the three predicted transcription factors (mef2, tcf1 and lef1) are linked to Wnt signaling pathway. It is already known that mef2 can enhance canonical Wnt signal [69] and it is involved in osteogenesis as well [70,71,72]. The binding site motif for mef2 appeared to be deleted in 3’UTR of the LT and riverine species. On the other hand, a binding site motif for tcf1 was gained in 3’UTR of the LM, LV and riverine species. In mice, tcf1 is demonstrated to be involve in paraxial mesoderm and limb formation and appeared to be act downstream of Wnt signal similar to lef1 transcription factor [73]. Moreover, canonical Wnt signaling has been shown to regulate osteogenesis through tcf1 responsive element on regulatory sequence of runx2 in mammals [74]. The third motif predicted to be a binding site for lef1 and only found within eda—tnfsf13b inter-genic region of LT species. Lef1 is again a well-known mediator of canonical Wnt signaling pathway which inhibits final stage of osteoblast differentiation [75] but it is essential for osteoblast proliferation and normal skeletal development [76, 77]. During development lef1 function is shown to be essential for scale outgrowth in zebrafish [78], and eda expression is known to be regulated by Wnt signal through lef1 transcriptional activity in mammals [79, 80]. In stickleback, mutation in an inter-genic region between eda and tnfsf13b genes is suggested to affect a binding site for c-jun transcription factor which its interaction with lef1 is required for eda transcriptional response to Wnt signal during armor plate formation [68]. Taken together, these observations, suggest mutations in enhancer sequences required for binding of Wnt signal components as potential underlying reason for the divergent expression of eda in both anterior and posterior scales of the cichlid species in this study.

In the posterior scale, in addition to eda, three more genes, dlx5, rankl and shh, displayed expression correlation with morphological divergence across the cichlid species (Fig. 6). The first gene, distal less homeobox 5 or dlx5, encodes transcription factor stimulating osteoblast differentiation and bone development, and it is also implicated in scale development and regeneration in fish [81,82,83]. Apart from its role in skeletogenesis, dlx5 has been found to be involved in divergent development and morphogenesis of other tissues in cichlids such as teeth and nuchal hump [84,85,86]. In goldfish, dlx5 expression appeared to be important at early stages of scale regeneration [81], and in both zebrafish and goldfish, dlx5 transcription in scale can be affected by environmental clues such as mechanical stimulus [82, 83]. The second gene, rankl, encodes a ligand for osteoprotegerin (opg) and plays a crucial role in osteoclast differentiation and bone remodelling. Changes in rankl transcription appeared to be important during scale regeneration in goldfish [81, 87], as well as intercellular communications regulating scale bone remodeling in zebrafish and goldfish [82, 88]. Both dlx5 and rankl have shown expression correlation patterns opposite to eda and shh in the posterior scales. Although, direct regulatory connections between these factors have not been investigated in scale but these findings suggest their potential interactions at transcriptional level. Moreover, higher expression of rankl in the scales of LT species compared to LM and LV species might indicate higher level of bone remodelling in their scales.

The third gene, sonic hedgehog or shh, encodes a ligand of hedgehog signaling pathway which is shown to control scale morphogenesis in relationship with the formation of the epidermal fold in the posterior region of scale in fish [16]. In zebrafish, epidermal expression of shh has been shown to regulate scale regeneration through controlling osteoblast population and affecting directional bone growth [89]. We found similar expression pattern between eda and shh which is more pronounced in the posterior scales. This is consistent with previous findings in other vertebrates, for instance, eda has been demonstrated to act upstream of shh and induce shh expression during ectodermal organogenesis in mammals (e.g. during hair placode formation) [90,91,92,93,94]. Furthermore, it has been shown that the eda-dependent regulation of shh might be a part of larger molecular cascade in which an upstream signal such as Wnt pathway activates eda signal and in turn eda induces shh transcription [92, 94, 95]. These observations suggest potential role of Wnt-eda-shh axis in divergent scale morphogenesis across Haplochromine cichlids, which seems to be more pronounced in the posterior scales. Finally, it is important to emphasize that since the expression differences are only investigated in mature stage, well beyond the larval stages at which scale formation initiates, it is likely that these findings mainly explain differences in growth process than morphogenesis of the scales.

Conclusions

This is the first attempt to study cross-species association between gene expression and morphological divergence in scales of cichlids from different lakes. Our results provide evidence for potential role of a key signal mediated by eda gene to be involved in divergent morphogenesis of scale in closely related cichlid species. We show that eda expression has lower level in the scales of species from the older lake (Lake Tanganyika) and correlates with the observed shape variations across species. Our findings shed light on molecular basis of morphological divergence of a less studied skeletal element; however, further transcriptional and functional investigations during scale development and growth/mineralization are required to assure eda signaling underlies the identified morphological differences and whether these differences have adaptive relevance in ecological and evolutionary-developmental contexts.

Availability of data and materials

All the gene expression data generated or analysed during this study are included in this published article. The transcriptome sequences used to design qPCR primers and check the variations in 3’UTR of eda genes for this article are available in the Sequencing Read Archive (SRA) of NCBI at https://www.ncbi.nlm.nih.gov/ and can be accessed through; A. burtoni (SRX4523155), L. trewavasae (SRX6432658), M. zebra (SRX8892877), N. omnicaeruleus (SRX8567938), P. sauvagei (SRX8892920), P. famula (SRX6445829), P. polyodon (ERX3501455), S. benthicola (ERX1818621), S. diagramma (SRX8567938) and T. tropheops (ERX659441). The genomic sequences used to check the variations in the inter-genic region of eda genes for this article are available in European Nucleotide Archive (ENA) at https://www.ebi.ac.uk/ena/ and can be accessed through; A. burtoni (GCA_000239415.1), P. famula (GCA_015108095.1), N. omnicaeruleus (SRR12700904), P. polyodon (GCA_015103895.1), S. fryeri (ERX1818621), S. diagramma (GCA_900408965.1), M. zebra (GCA_000238955.1), T. tropheops (SAMEA2661272), L. trewavasae (SAMN12216683), and P. sauvagei (GCA_018403495.1).

Abbreviations

LT:

Lake Tanganyika

LM:

Lake Malawi

LV:

Lake Victoria

bmp4 :

Bone morphogenetic protein 4

col1a2 :

Collagen type I alpha 2 chain

ctsk :

Cathepsin K

dlx5 :

Distal-less homeobox 5

eda :

Ectodysplasin A

edar :

Ectodysplasin A receptor

fgf20 :

Fibroblast growth factor 20

fgfr1 :

Fibroblast growth factor receptor 1

mmp2 :

Matrix metallopeptidase 2

mmp9 :

Matrix metallopeptidase 9

opg :

Osteoprotegerin

rankl :

Receptor activator of nuclear factor kappa B ligand

runx2a :

Runt-related transcription factor 2 alpha

sema4d :

Semaphorin 4D

shh :

Sonic hedgehog signaling molecule

sp7 :

Osterix transcription factor

References

  1. Fryer G, Iles TD. The cichlid fishes of the great lakes of Africa: their biology and evolution. Edinburgh: Oliver and Boyd; 1972.

    Google Scholar 

  2. Kocher TD. Adaptive evolution and explosive speciation: the cichlid fish model. Nat Rev Genet. 2004. https://doi.org/10.1038/nrg1316.

    Article  PubMed  Google Scholar 

  3. Salzburger W. Understanding explosive diversification through cichlid fish genomics. Nat Rev Genet. 2018. https://doi.org/10.1038/s41576-018-0043-9.

    Article  PubMed  Google Scholar 

  4. Muschick M, Indermaur A, Salzburger W. Convergent evolution within an adaptive radiation of cichlid fishes. Curr Biol. 2012;22:2362–8.

    Article  CAS  PubMed  Google Scholar 

  5. Wanek AK, Sturmbauer C. Form, function and phylogeny: comparative morphometrics of Lake Tanganyika’s cichlid tribe Tropheini. Zool Scr. 2015;44:362–73.

    Article  PubMed  PubMed Central  Google Scholar 

  6. Albertson RC, Kocher TD. Genetic and developmental basis of cichlid trophic diversity. Heredity (Edinb). 2006;97:211–21. https://doi.org/10.1038/sj.hdy.6800864.

    Article  CAS  Google Scholar 

  7. Singh P, Börger C, More H, Sturmbauer C. The role of alternative splicing and differential gene expression in cichlid adaptive radiation. Genome Biol Evol. 2017;9:2764–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. Albertson RC, Streelman JT, Kocher TD, Yelick PC. Integration and evolution of the cichlid mandible: the molecular basis of alternate feeding strategies. Proc Natl Acad Sci U S A. 2005;102:16287–92. https://doi.org/10.1073/pnas.0506649102.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Muschick M, Barluenga M, Salzburger W, Meyer A. Adaptive phenotypic plasticity in the Midas cichlid fish pharyngeal jaw and its relevance in adaptive radiation. BMC Evol Biol. 2011;11:116. https://doi.org/10.1186/1471-2148-11-116.

    Article  PubMed  PubMed Central  Google Scholar 

  10. Ahi EP, Singh P, Duenser A, Gessl W, Sturmbauer C. Divergence in larval jaw gene expression reflects differential trophic adaptation in haplochromine cichlids prior to foraging. BMC Evol Biol. 2019;19:150. https://doi.org/10.1186/s12862-019-1483-3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Singh P, Ahi EP, Sturmbauer C. Gene coexpression networks reveal molecular interactions underlying cichlid jaw modularity. BMC Ecol Evol. 2021;21:1–17. https://doi.org/10.1186/s12862-021-01787-9.

    Article  Google Scholar 

  12. Parsons KJ, Concannon M, Navon D, Wang J, Ea I, Groveas K, et al. Foraging environment determines the genetic architecture and evolutionary potential of trophic morphology in cichlid fishes. Mol Ecol. 2016. https://doi.org/10.1111/mec.13801.

    Article  PubMed  Google Scholar 

  13. van Rijssel JC, Hoogwater ES, Kishe-Machumu MA, van Reenen E, Spits KV, van der Stelt RC, et al. Fast adaptive responses in the oral jaw of Lake Victoria cichlids. Evolution (N Y). 2015;69:179–89. https://doi.org/10.1111/evo.12561.

    Article  Google Scholar 

  14. Huysseune A, Sire JY. Evolution of patterns and processes in teeth and tooth-related tissues in non-mammalian vertebrates. Eur J Oral Sci. 1998. https://doi.org/10.1111/j.1600-0722.1998.tb02211.x.

    Article  PubMed  Google Scholar 

  15. Sire JY, Huysseune A. Formation of dermal skeletal and dental tissues in fish: A comparative and evolutionary approach. Biol Rev. 2003. https://doi.org/10.1017/S1464793102006073.

    Article  PubMed  Google Scholar 

  16. Sire JY, Akimenko MA. Scale development in fish: a review, with description of sonic hedgehog (shh) expression in the zebrafish (Danio rerio). Int J Dev Biol. 2004;48:233–47.

    Article  CAS  PubMed  Google Scholar 

  17. Wainwright DK, Lauder G V. Mucus Matters: The Slippery and Complex Surfaces of Fish. 2017. https://doi.org/10.1007/978-3-319-74144-4_10

  18. Kondo S, Kuwahara Y, Kondo M, Naruse K, Mitani H, Wakamatsu Y, et al. The medaka rs-3 locus required for scale development encodes ectodysplasin-A receptor. Curr Biol. 2001;11:1202–6.

    Article  CAS  PubMed  Google Scholar 

  19. Harris MP, Rohner N, Schwarz H, Perathoner S, Konstantinidis P, Nüsslein-Volhard C. Zebrafish eda and edar Mutants Reveal Conserved and Ancestral Roles of Ectodysplasin Signaling in Vertebrates. PLoS Genet. 2008;4: e1000206. https://doi.org/10.1371/journal.pgen.1000206.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Albertson RC, Kawasaki KC, Tetrault ER, Powder KE. Genetic analyses in Lake Malawi cichlids identify new roles for Fgf signaling in scale shape variation. Commun Biol. 2018;1:1–11.

    Article  Google Scholar 

  21. Biggs LC, Mikkola ML. Early inductive events in ectodermal appendage morphogenesis. Seminars Cell Dev Biol. 2014. https://doi.org/10.1016/j.semcdb.2014.01.007.

    Article  Google Scholar 

  22. Mikkola ML, Thesleff I. Ectodysplasin signaling in development. Cytokine Growth Factor Rev. 2003;14:211–24.

    Article  CAS  PubMed  Google Scholar 

  23. Colosimo PF, Hosemann KE, Balabhadra S, Villarreal G, Dickson H, Grimwood J, et al. Widespread parallel evolution in sticklebacks by repeated fixation of ectodysplasin alleles. Science. 2005. https://doi.org/10.1126/science.1107239.

    Article  PubMed  Google Scholar 

  24. Viertler A, Salzburger W, Ronco F. Comparative scale morphology in the adaptive radiation of cichlid fishes (Perciformes: Cichlidae) from Lake Tanganyika. Biol J Linn Soc. 2021. https://doi.org/10.1093/biolinnean/blab099.

    Article  Google Scholar 

  25. Rohlf FJ, Slice D. Extensions of the procrustes method for the optimal superimposition of landmarks. Syst Zool. 1990;39:40. https://doi.org/10.2307/2992207.

    Article  Google Scholar 

  26. Klingenberg CP. MorphoJ: an integrated software package for geometric morphometrics. Mol Ecol Resour. 2011;11:353–7. https://doi.org/10.1111/j.1755-0998.2010.02924.x.

    Article  PubMed  Google Scholar 

  27. Hammer Ø, Harper DAT, Ryan PD. Past: Paleontological statistics software package for education and data analysis. Palaeontol Electron. 2001;4(1):1–9.

    Google Scholar 

  28. Mitteroecker P, Gunz P. Advances in Geometric morphometrics. Evol Biol. 2009;36:235–47.

    Article  Google Scholar 

  29. Team RC. A language and environment for statistical computing. Vienna, Austria: R Foundation for Statistical Computing; 2018. http://www.R-project.org.

  30. Yang CG, Wang XL, Tian J, Liu W, Wu F, Jiang M, et al. Evaluation of reference genes for quantitative real-time RT-PCR analysis of gene expression in Nile tilapia (Oreochromis niloticus). Gene. 2013;527:183–92.

    Article  CAS  PubMed  Google Scholar 

  31. Gunter HM, Fan S, Xiong F, Franchini P, Fruciano C, Meyer A. Shaping development through mechanical strain: the transcriptional basis of diet-induced phenotypic plasticity in a cichlid fish. Mol Ecol. 2013;22:4516–31.

    Article  PubMed  Google Scholar 

  32. Ahi EP, Richter F, Sefc KM. A gene expression study of ornamental fin shape in Neolamprologus brichardi, an African cichlid species. Sci Rep. 2017. https://doi.org/10.1038/s41598-017-17778-0.

    Article  PubMed  PubMed Central  Google Scholar 

  33. Ahi EP, Duenser A, Singh P, Gessl W, Sturmbauer C. Appetite regulating genes may contribute to herbivory versus carnivory trophic divergence in haplochromine cichlids. PeerJ. 2020;8: e8375. https://doi.org/10.7717/peerj.8375.

    Article  PubMed  PubMed Central  Google Scholar 

  34. Ahi EP, Richter F, Lecaudey LA, Sefc KM. Gene expression profiling suggests differences in molecular mechanisms of fin elongation between cichlid species. Sci Rep. 2019. https://doi.org/10.1038/s41598-019-45599-w.

    Article  PubMed  PubMed Central  Google Scholar 

  35. Lecaudey LA, Sturmbauer C, Singh P, Ahi EP. Molecular mechanisms underlying nuchal hump formation in dolphin cichlid. Cyrtocara moorii Sci Rep. 2019;9:20296. https://doi.org/10.1038/s41598-019-56771-7.

    Article  CAS  PubMed  Google Scholar 

  36. Brawand D, Wagner CE, Li YI, Malinsky M, Keller I, Fan S, et al. The genomic substrate for adaptive radiation in African cichlid fish. Nature. 2014;513:375–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Lecaudey LA, Singh P, Sturmbauer C, Duenser A, Gessl W, Ahi EP. Transcriptomics unravels molecular players shaping dorsal lip hypertrophy in the vacuum cleaner cichlid. Gnathochromis permaxillaris BMC Genomics. 2021;22:506.

    Article  CAS  PubMed  Google Scholar 

  38. Ahi EP, Lecaudey LA, Ziegelbecker A, Steiner O, Glabonjat R, Goessler W, et al. Comparative transcriptomics reveals candidate carotenoid color genes in an East African cichlid fish. BMC Genomics. 2020. https://doi.org/10.1186/s12864-020-6473-8.

    Article  PubMed  PubMed Central  Google Scholar 

  39. Ahi EP, Lecaudey LA, Ziegelbecker A, Steiner O, Goessler W, Sefc KM. Expression levels of the tetratricopeptide repeat protein gene ttc39b covary with carotenoid-based skin colour in cichlid fish. Biol Lett. 2020. https://doi.org/10.1098/rsbl.2020.0629.

    Article  PubMed  PubMed Central  Google Scholar 

  40. Zerbino DR, Achuthan P, Akanni W, Amode MR, Barrell D, Bhai J, et al. Ensembl 2018. Nucleic Acids Res. 2018;46:D754–61.

    Article  CAS  PubMed  Google Scholar 

  41. Fleige S, Pfaffl MW. RNA integrity and the effect on the real-time qRT-PCR performance. Mol Aspects Med. 2006;27:126–39.

    Article  CAS  PubMed  Google Scholar 

  42. Mahony S, Benos PV. STAMP: a web tool for exploring DNA-binding motif similarities. Nucleic Acids Res. 2007. https://doi.org/10.1093/nar/gkm272.

    Article  PubMed  PubMed Central  Google Scholar 

  43. Matys V, Fricke E, Geffers R, Gößling E, Haubrock M, Hehl R, et al. TRANSFAC®: Transcriptional regulation, from patterns to profiles. Nucleic Acids Res. 2003;31:374–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  44. Hellemans J, Mortier G, De Paepe A, Speleman F, Vandesompele J. qBase relative quantification framework and software for management and automated analysis of real-time quantitative PCR data. Genome Biol. 2007;8:R19. https://doi.org/10.1186/gb-2007-8-2-r19.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Pashay Ahi E, Sefc KM. Towards a gene regulatory network shaping the fins of the Princess cichlid. Sci Rep. 2018;8:9602. https://doi.org/10.1038/s41598-018-27977-y.

    Article  CAS  Google Scholar 

  46. Pfaffl MW, Tichopad A, Prgomet C, Neuvians TP. Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper–Excel-based tool using pair-wise correlations. Biotechnol Lett. 2004;26:509–15.

    Article  CAS  PubMed  Google Scholar 

  47. Andersen CL, Jensen JL, Ørntoft TF. Normalization of real-time quantitative reverse transcription-PCR data: a model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets. Cancer Res. 2004;64:5245–50. https://doi.org/10.1158/0008-5472.CAN-04-0496.

    Article  CAS  PubMed  Google Scholar 

  48. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, et al. Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biol. 2002. https://doi.org/10.1186/gb-2002-3-7-research0034.

    Article  PubMed  PubMed Central  Google Scholar 

  49. Pfaffl MW. A new mathematical model for relative quantification in real-time RT-PCR. Nucleic Acids Res. 2001;29:e45.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  50. Bergkvist A, Rusnakova V, Sindelka R, Garda JMA, Sjögreen B, Lindh D, et al. Gene expression profiling—clusters of possibilities. Methods. 2010;50:323–35. https://doi.org/10.1016/j.ymeth.2010.01.009.

    Article  CAS  PubMed  Google Scholar 

  51. Kubista M, Andrade JM, Bengtsson M, Forootan A, Jonák J, Lind K, et al. The real-time polymerase chain reaction. Mol Aspects Med. 2006;27:95–125. https://doi.org/10.1016/j.mam.2005.12.007.

    Article  CAS  PubMed  Google Scholar 

  52. Ronco F, Matschiner M, Böhne A, Boila A, Büscher HH, El Taher A, et al. Drivers and dynamics of a massive adaptive radiation in cichlid fishes. Nature. 2021;589:76–81.

    Article  CAS  PubMed  Google Scholar 

  53. Salzburger W, Meyer A, Baric S, Verheyen E, Sturmbauer C. Phylogeny of the Lake Tanganyika cichlid species flock and its relationship to the Central and East African haplochromine cichlid fish faunas. Syst Biol. 2002. https://doi.org/10.1080/106351502753475907.

    Article  PubMed  Google Scholar 

  54. Koblmüller S, Schliewen UK, Duftner N, Sefc KM, Katongo C, Sturmbauer C. Age and spread of the haplochromine cichlid fishes in Africa. Mol Phylogenet Evol. 2008;49:153–69. https://doi.org/10.1016/J.YMPEV.2008.05.045.

    Article  PubMed  Google Scholar 

  55. Cohen AS, Soreghan MJ, Scholz CA. Estimating the age of formation of lakes: an example from Lake Tanganyika. East African Rift system Geology. 1993;21:511–4.

    CAS  Google Scholar 

  56. Delvaux D. Age of Lake Malawi (Nyasa) and water level fluctuations. Mus roy Afrr centr, Tervuren (Belg), Dept Geol Min, Rapp ann 1993 1994. 1995;108.

  57. Johnson TC, Kelts K, Odada E. The holocene history of lake victoria. Ambio. 2000. https://doi.org/10.1579/0044-7447-29.1.2.

    Article  Google Scholar 

  58. Blažek R, Polačik M, Smith C, Honza M, Meyer A, Reichard M. Success of cuckoo catfish brood parasitism reflects coevolutionary history and individual experience of their cichlid hosts. Sci Adv. 2018. https://doi.org/10.1126/sciadv.aar4380.

    Article  PubMed  PubMed Central  Google Scholar 

  59. Koblmüller S, Duftner N, Sefc KM, Aibara M, Stipacek M, Blanc M, et al. Reticulate phylogeny of gastropod-shell-breeding cichlids from Lake Tanganyika—the result of repeated introgressive hybridization. BMC Evol Biol. 2007. https://doi.org/10.1186/1471-2148-7-7.

    Article  PubMed  PubMed Central  Google Scholar 

  60. Takamura K. Interspecific relationships of aufwuchs-eating fishes in Lake Tanganyika. Environ Biol Fishes. 1984;10:225–41.

    Article  Google Scholar 

  61. Takahashi R, Watanabe K, Nishida M, Hori M. Evolution of feeding specialization in Tanganyikan scale-eating cichlids: a molecular phylogenetic approach. BMC Evol Biol. 2007. https://doi.org/10.1186/1471-2148-7-195.

    Article  PubMed  PubMed Central  Google Scholar 

  62. Kovac R, Boileau N, Muschick M, Salzburger W. The diverse prey spectrum of the Tanganyikan scale-eater Perissodus microlepis (Boulenger, 1898). Hydrobiologia. 2019. https://doi.org/10.1007/s10750-018-3714-9.

    Article  Google Scholar 

  63. Boileau N, Cortesi F, Egger B, Muschick M, Indermaur A, Theis A, et al. A complex mode of aggressive mimicry in a scale-eating cichlid fish. Biol Lett. 2015. https://doi.org/10.1098/rsbl.2015.0521.

    Article  PubMed  PubMed Central  Google Scholar 

  64. Irisarri I, Singh P, Koblmüller S, Torres-Dowdall J, Henning F, Franchini P, et al. Phylogenomics uncovers early hybridization and adaptive loci shaping the radiation of Lake Tanganyika cichlid fishes. Nat Commun. 2018;9:1–12.

    Article  CAS  Google Scholar 

  65. Mikkola ML. Molecular aspects of hypohidrotic ectodermal dysplasia. Am J Med Genet Part A. 2009;149A:2031–6. https://doi.org/10.1002/ajmg.a.32855.

    Article  CAS  PubMed  Google Scholar 

  66. Cheng J, Sedlazek F, Altmüller J, Nolte AW. Ectodysplasin signalling genes and phenotypic evolution in sculpins ( Cottus ). Proc R Soc B Biol Sci. 2015;282:20150746. https://doi.org/10.1098/rspb.2015.0746.

    Article  CAS  Google Scholar 

  67. Shono T, Thiery AP, Cooper RL, Kurokawa D, Britz R, Okabe M, et al. Evolution and developmental diversity of skin spines in pufferfishes. iScience. 2019;19:1248–59.

    Article  PubMed  PubMed Central  Google Scholar 

  68. O’brown NM, Summers BR, Jones FC, Brady SD, Kingsley DM. A recurrent regulatory change underlying altered expression and Wnt response of the stickleback armor plates gene EDA. Elife. 2015. https://doi.org/10.7554/eLife.05290.

    Article  PubMed  PubMed Central  Google Scholar 

  69. Ehyai S, Dionyssiou MG, Gordon JW, Williams D, Siu KWM, McDermott JC. A p38 MAPK regulated MEF2:β-catenin interaction enhances canonical Wnt signalling. Mol Cell Biol. 2015. https://doi.org/10.1128/MCB.00832-15.

    Article  PubMed  Google Scholar 

  70. Shen S, Huang D, Feng G, Zhu L, Zhang Y, Cao P, et al. MEF2 transcription factor regulates osteogenic differentiation of dental pulp stem cells. Cell Reprogram. 2016;18:237–45. https://doi.org/10.1089/cell.2016.0016.

    Article  CAS  PubMed  Google Scholar 

  71. Kawane T, Komori H, Liu W, Moriishi T, Miyazaki T, Mori M, et al. Dlx5 and Mef2 regulate a novel Runx2 enhancer for osteoblast-specific expression. J Bone Miner Res. 2014;29:1960–9. https://doi.org/10.1002/jbmr.2240.

    Article  CAS  PubMed  Google Scholar 

  72. Potthoff MJ, Olson EN. MEF2: a central regulator of diverse developmental programs. Development. 2007;134:4131–40.

    Article  CAS  PubMed  Google Scholar 

  73. Galceran J, Fariñas I, Depew MJ, Clevers H, Grosschedl R. Wnt3a(−/−) -like phenotype and limb deficiency in Lef1(−/−)Tcf1(−/−) mice. Genes Dev. 1999;13:709–17.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Gaur T, Lengner CJ, Hovhannisyan H, Bhat RA, Bodine PVN, Komm BS, et al. Canonical WNT signaling promotes osteogenesis by directly stimulating Runx2 gene expression. J Biol Chem. 2005;280:33132–40. https://doi.org/10.1074/jbc.M500608200.

    Article  CAS  PubMed  Google Scholar 

  75. Kahler RA, Galindo M, Lian J, Stein GS, van Wijnen AJ, Westendorf JJ. Lymphocyte enhancer-binding factor 1 (Lef1) inhibits terminal differentiation of osteoblasts. J Cell Biochem. 2006;97:969–83. https://doi.org/10.1002/jcb.20702.

    Article  CAS  PubMed  Google Scholar 

  76. Kahler RA, Westendorf JJ. Lymphoid enhancer factor-1 and β-catenin inhibit Runx2-dependent transcriptional activation of the osteocalcin promoter. J Biol Chem. 2003;278:11937–44.

    Article  CAS  PubMed  Google Scholar 

  77. Hoeppner LH, Secreto F, Jensen ED, Li X, Kahler RA, Westendorf JJ. Runx2 and bone morphogenic protein 2 regulate the expression of an alternative Lef1 transcript during osteoblast maturation. J Cell Physiol. 2009;221:480–9. https://doi.org/10.1002/jcp.21879.

    Article  CAS  PubMed  Google Scholar 

  78. Aman AJ, Fulbright AN, Parichy DM. Wnt/β-catenin regulates an ancient signaling network during zebrafish scale development. Elife. 2018. https://doi.org/10.7554/eLife.37001.

    Article  PubMed  PubMed Central  Google Scholar 

  79. Durmowicz MC, Cui CY, Schlessinger D. The EDA gene is a target of, but does not regulate Wnt signaling. Gene. 2002;285:203–11.

    Article  CAS  PubMed  Google Scholar 

  80. Lévy J, Capri Y, Rachid M, Dupont C, Vermeesch JR, Devriendt K, et al. LEF1 haploinsufficiency causes ectodermal dysplasia. Clin Genet. 2020;97:595–600. https://doi.org/10.1111/cge.13714.

    Article  CAS  PubMed  Google Scholar 

  81. Thamamongood TA, Furuya R, Fukuba S, Nakamura M, Suzuki N, Hattori A. Expression of osteoblastic and osteoclastic genes during spontaneous regeneration and autotransplantation of goldfish scale: a new tool to study intramembranous bone regeneration. Bone. 2012;50:1240–9.

    Article  CAS  PubMed  Google Scholar 

  82. Kitamura K, Takahira K, Inari M, Satoh Y, Hayakawa K, Tabuchi Y, et al. Zebrafish scales respond differently to in vitro dynamic and static acceleration: analysis of interaction between osteoblasts and osteoclasts. Comp Biochem Physiol A Mol Integr Physiol. 2013;166:74–80.

    Article  CAS  PubMed  Google Scholar 

  83. Suzuki N, Hanmoto T, Yano S, Furusawa Y, Ikegame M, Tabuchi Y, et al. Low-intensity pulsed ultrasound induces apoptosis in osteoclasts: fish scales are a suitable model for the analysis of bone metabolism by ultrasound. Comp Biochem Physiol -Part A Mol Integr Physiol. 2016;195:26–31.

    Article  CAS  Google Scholar 

  84. Renz AJ, Gunter HM, Fischer JMF, Qiu H, Meyer A, Kuraku S. Ancestral and derived attributes of the dlx gene repertoire, cluster structure and expression patterns in an African cichlid fish. EvoDevo. 2011;2:1. https://doi.org/10.1186/2041-9139-2-1.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  85. Hulsey CD, Fraser GJ, Meyer A. Biting into the genome to phenome map: developmental genetic modularity of cichlid fish dentitions. Integr Comp Biol. 2016. https://doi.org/10.1093/icb/icw059.

    Article  PubMed  Google Scholar 

  86. Lecaudey LA, Sturmbauer C, Singh P, Ahi EP. Molecular mechanisms underlying nuchal hump formation in dolphin cichlid Cyrtocara moorii. Sci Rep. 2019;9:20296.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  87. Yamamoto T, Ikegame M, Kawago U, Tabuchi Y, Hirayama J, Sekiguchi T, et al. Detection of RANKL-producing cells and osteoclastic activation by the addition of exogenous RANKL in the regenerating scales of goldfish. Biol Sci Sp. 2020;34:34–40.

    Article  Google Scholar 

  88. Tazaki Y, Sugitani K, Ogai K, Kobayashi I, Kawasaki H, Aoyama T, et al. RANKL, Ephrin-Eph and Wnt10b are key intercellular communication molecules regulating bone remodeling in autologous transplanted goldfish scales. Comp Biochem Physiol Part A Mol Integr Physiol. 2018;225:46–58.

    Article  CAS  Google Scholar 

  89. Iwasaki M, Kuroda J, Kawakami K, Wada H. Epidermal regulation of bone morphogenesis through the development and regeneration of osteoblasts in the zebrafish scale. Dev Biol. 2018;437:105–19.

    Article  CAS  PubMed  Google Scholar 

  90. Schmidt-Ullrich R, Tobin DJ, Lenhard D, Schneider P, Paus R, Scheidereit C. NF-κB transmits Eda A1/EdaR signalling to activate Shh and cyclin D1 expression, and controls post-initiation hair placode down growth. Development. 2006;133:1045–57.

    Article  CAS  PubMed  Google Scholar 

  91. Fliniaux I, Mikkola ML, Lefebvre S, Thesleff I. Identification of dkk4 as a target of Eda-A1/Edar pathway reveals an unexpected role of ectodysplasin as inhibitor of Wnt signalling in ectodermal placodes. Dev Biol. 2008;320:60–71.

    Article  CAS  PubMed  Google Scholar 

  92. Cui CY, Yin M, Sima J, Childress V, Michel M, Piao Y, et al. Involvement of Wnt, Eda and Shh at defined stages of sweat gland development. Dev. 2014;141:3752–60.

    Article  CAS  Google Scholar 

  93. Pummila M, Fliniaux I, Jaatinen R, James MJ, Laurikkala J, Schneider P, et al. Ectodysplasin has a dual role in ectodermal organogenesis: inhibition of Bmp activity and induction of Shh expression. Development. 2007;134:117–25.

    Article  CAS  PubMed  Google Scholar 

  94. Xiao Y, Thoresen DT, Miao L, Williams JS, Wang C, Atit RP, et al. A cascade of Wnt, Eda, and Shh signaling is essential for touch dome merkel cell development. PLOS Genet. 2016;12: e1006150. https://doi.org/10.1371/journal.pgen.1006150.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  95. Hammerschmidt B, Schlake T. Localization of Shh expression by Wnt and Eda affects axial polarity and shape of hairs. Dev Biol. 2007;305:246–61.

    Article  CAS  PubMed  Google Scholar 

  96. Hong YJ, Choi YW, Myung KB, Choi HY. The immunohistochemical patterns of calcification-related molecules in the epidermis and dermis of the zebrafish (Danio rerio). Ann Dermatol. 2011;23:299–303.

    Article  PubMed  PubMed Central  Google Scholar 

  97. Le Guellec D, Morvan-Dubois G, Sire JY. Skin development in bony fish with particular emphasis on collagen deposition in the dermis of the zebrafish (Danio rerio). Int J Dev Biol. 2004;48:217–31.

    Article  PubMed  Google Scholar 

  98. Metz JR, Leeuwis RHJ, Zethof J, Flik G. Zebrafish ( Danio rerio ) in calcium-poor water mobilise calcium and phosphorus from scales. J Appl Ichthyol. 2014;30:671–7.

    Article  CAS  Google Scholar 

  99. de Vrieze E, Metz JR, Von den Hoff JW, Flik G. ALP, TRAcP and cathepsin K in elasmoid scales: a role in mineral metabolism? J Appl Ichthyol. 2010;26:210–3.

    Article  Google Scholar 

  100. Iida Y, Hibiya K, Inohaya K, Kudo A. Eda/Edar signaling guides fin ray formation with preceding osteoblast differentiation, as revealed by analyses of the medaka all-fin less mutant afl. Dev Dyn. 2014;243:765–77. https://doi.org/10.1002/dvdy.24120.

    Article  CAS  PubMed  Google Scholar 

  101. Atukorala ADS, Inohaya K, Baba O, Tabata MJ, Ratnayake RAR, Abduweli D, et al. Scale and tooth phenotypes in medaka with a mutated ectodysplasin-A receptor: implications for the evolutionary origin of oral and pharyngeal teeth. Arch Histol Cytol. 2011;73:139–48. https://doi.org/10.1679/aohc.73.139.

    Article  Google Scholar 

  102. Daane JM, Rohner N, Konstantinidis P, Djuranovic S, Harris MP. Parallelism and epistasis in skeletal evolution identified through use of phylogenomic mapping strategies. Mol Biol Evol. 2016;33:162–73. https://doi.org/10.1093/molbev/msv208.

    Article  CAS  PubMed  Google Scholar 

  103. Cooper WJ, Wirgau RM, Sweet EM, Albertson RC. Deficiency of zebrafish fgf20a results in aberrant skull remodeling that mimics both human cranial disease and evolutionarily important fish skull morphologies. Evol Dev. 2013;15:426–41. https://doi.org/10.1111/ede.12052.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  104. Rohner N, Bercsényi M, Orbán L, Kolanczyk ME, Linke D, Brand M, et al. Duplication of fgfr1 Permits Fgf Signaling to Serve as a Target for Selection during Domestication. Curr Biol. 2009;19:1642–7.

    Article  CAS  PubMed  Google Scholar 

  105. de Vrieze E, Sharif F, Metz JR, Flik G, Richardson MK. Matrix metalloproteinases in osteoclasts of ontogenetic and regenerating zebrafish scales. Bone. 2011;48:704–12.

    Article  PubMed  Google Scholar 

  106. Iimura K, Tohse H, Ura K, Takagi Y. Expression patterns of runx2, sparc, and bgp during scale regeneration in the goldfish carassius auratus. J Exp Zool Part B Mol Dev Evol. 2012;318:190–8.

    Article  CAS  Google Scholar 

  107. Zhong Z, Niu P, Wang M, Huang G, Xu S, Sun Y, et al. Targeted disruption of sp7 and myostatin with CRISPR-Cas9 results in severe bone defects and more muscular cells in common carp. Sci Rep. 2016;6:1–14.

    Google Scholar 

Download references

Acknowledgements

The authors acknowledge Institute of Biology at University of Graz for providing fish breeding and laboratory facilities, and the Austrian Science Fund for the financial support. We also thank Stephan Koblmüller for sharing his precious knowledge on cichlid fishes of the African Lakes.

Funding

This study was funded by the Austrian Science Fund (Grant P29838). The Austrian Science Fund requires clarification of all legal issues concerning animal keeping, animal experiments and sampling design prior to grant submission and evaluation, but does not interfere in writing and data interpretation, but funds open access of the resulting publications.

Author information

Author notes

  1. Maximilian Wagner and Sandra Bračun contributed equally to this work

Authors and Affiliations

  1. Institute of Biology, University of Graz, Universitätsplatz 2, 8010, Graz, Austria

    Maximilian Wagner, Sandra Bračun, Anna Duenser, Christian Sturmbauer, Wolfgang Gessl & Ehsan Pashay Ahi

  2. Department of Biology, University of Antwerp, Groenenborgerlaan 171, 2020, Antwerp, Belgium

    Maximilian Wagner

  3. Organismal and Evolutionary Biology Research Programme, University of Helsinki, Viikinkaari 9, 00014, Helsinki, Finland

    Ehsan Pashay Ahi

Authors
  1. Maximilian Wagner
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Sandra Bračun
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Anna Duenser
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Christian Sturmbauer
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Wolfgang Gessl
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Ehsan Pashay Ahi
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

EPA, SB, MW and CS designed the study. SB, EPA, MW, and AD conducted the laboratory experiment, measurements and figure preparations. MW and EPA analysed the data, and EPA, MW, AD and CS wrote the manuscript. WG and AD performed fish breeding and sampling. WG photographed the adult fishes used in Fig. 1A. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Christian Sturmbauer or Ehsan Pashay Ahi.

Ethics declarations

Ethics approval and consent to participate

All experimental protocols were approved by the Federal Ministry of Science, Research and Economy of Austria (BMWFW) under permit BMWFW-66.007/0004-WF/V/3b/2016. All methods were carried out in accordance with relevant guidelines and regulations of the Austrian animal welfare law. All methods are reported in accordance with ARRIVE guidelines (https://arriveguidelines.org) for the reporting of animal experiments.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

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

Supplementary Information

Additional file 1.

Information about qPCR primers and raw expression data.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Wagner, M., Bračun, S., Duenser, A. et al. Expression variations in ectodysplasin-A gene (eda) may contribute to morphological divergence of scales in haplochromine cichlids. BMC Ecol Evo 22, 28 (2022). https://doi.org/10.1186/s12862-022-01984-0

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12862-022-01984-0

Keywords

BMC Ecology and Evolution

ISSN: 2730-7182

Contact us