A Genomic Cluster Containing Novel and Conserved Genes is Associated with Cichlid Fish Dental Developmental Convergence

Abstract

The two toothed jaws of cichlid fishes provide textbook examples of convergent evolution. Tooth phenotypes such as enlarged molar-like teeth used to process hard-shelled mollusks have evolved numerous times independently during cichlid diversification. Although the ecological benefit of molar-like teeth to crush prey is known, it is unclear whether the same molecular mechanisms underlie these convergent traits. To identify genes involved in the evolution and development of enlarged cichlid teeth, we performed RNA-seq on the serially homologous-toothed oral and pharyngeal jaws as well as the fourth toothless gill arch of Astatoreochromis alluaudi. We identified 27 genes that are highly upregulated on both tooth-bearing jaws compared with the toothless gill arch. Most of these genes have never been reported to play a role in tooth formation. Two of these genes (unk, rpfA) are not found in other vertebrate genomes but are present in all cichlid genomes. They also cluster genomically with two other highly expressed tooth genes (odam, scpp5) that exhibit conserved expression during vertebrate odontogenesis. Unk and rpfA were confirmed via in situ hybridization to be expressed in developing teeth of Astatotilapia burtoni. We then examined expression of the cluster’s four genes in six evolutionarily independent and phylogenetically disparate cichlid species pairs each with a large- and a small-toothed species. Odam and unk commonly and scpp5 and rpfA always showed higher expression in larger toothed cichlid jaws. Convergent trophic adaptations across cichlid diversity are associated with the repeated developmental deployment of this genomic cluster containing conserved and novel cichlid-specific genes.

Introduction

Convergent adaptations can evolve in response to similar environmental challenges (Losos 2011). However, whether these repeatedly evolved traits develop via novel or conserved genetic mechanisms is often unclear (Shapiro et al. 2004; Kratochwil et al. 2018). The evolution of similar traits could frequently arise through modifications of the expression of different genes (Härer et al. 2018). Alternatively, in some convergent traits such as the loss of stickleback pelvic spines (Shapiro et al. 2006) or acquired toxin resistance in newt-eating snakes (Geffeney et al. 2002), repeated modification of the same genetic loci are known to underlie convergent evolution. Attributes such as whether genes generating convergent traits are highly conserved (Joron et al. 2006), if gene expression shows limited pleiotropy (Papakostas et al. 2014), or how recently convergent species shared a common ancestor (Geffeney et al. 2002; Torres-Dowdall et al. 2017; Kratochwil et al. 2018) are all thought to contribute to whether similar versus novel developmental mechanisms underlie convergence. Yet, few systems with well-characterized genomes exhibit convergence in duplicated structures where shared expression is highly likely and also display evolutionary replication over different time-scales. However, the phenotypic convergence found in independent adaptive radiations of cichlid fishes provide natural evolutionary experiments that offer exceptional opportunities to study developmental genomic similarities in repeatedly evolved and ecologically relevant traits (Stiassny and Meyer 1999; Torres-Dowdall et al. 2017; Hulsey et al. 2018). For instance, phenotypes such as enlarged molar-like teeth that are used to process hard-shelled mollusks have evolved numerous times independently during cichlid diversification. Using an integrative combination of gene expression and comparative analyses, we examined the developmental genetic role of novel and conserved mechanisms during the convergent evolution of molariform teeth across cichlid fish diversity.

Teeth are thought to have evolved once in the ancestor of all jawed vertebrates (Smith and Coates 1998) making teeth a useful model for comparisons of vertebrate development. We also have extensive knowledge concerning the developmental genetic basis of odontogenesis from mammalian model systems of tooth development (Thesleff and Hurmerinta 1981; Chen et al. 1996; Helsinki 1996; Bei and Maas 1998; Dassule et al. 2000; Jernvall et al. 2000; Thesleff et al. 2001; Tucker and Sharpe 2004; Thesleff 2006; Wang et al. 2007; Popa et al. 2019). Importantly, many genes contributing to odontogenesis in mammals also play a role in tooth formation in other vertebrates as evolutionarily divergent as teleost fishes (Wise and Stock 2006; Harris et al. 2008; Vonk et al. 2008; Fraser et al. 2009, 2013; Ellis et al. 2015; Hulsey et al. 2016). But, since mammalian dentition differs substantially in phenotype from those of most other vertebrates, there could also be many novel mechanisms involved in tooth formation that were either lost in mammals or novelly emerged in teleosts. For example, mammals replace their teeth at most a single time (van Nievelt and Smith 2005; Luo 2007). However, cichlids and most vertebrates constantly replace their teeth throughout their lifetimes (Tuisku and Hildebrand 1994; Richman and Handrigan 2011; Fraser et al. 2013; Rasch et al. 2016). Furthermore, processes such as the teleost-specific genome duplication or gene birth from noncoding elements could have provided modern teleosts with novel odontogenetic loci (Ohno 1970; Kawasaki and Weiss 2003; Wu et al. 2011). Thus, tooth gene diversity in cichlids and other teleosts is likely to consist of both a core set of genes shared across most vertebrates as well as a potentially diverse group of genes specific to particular lineages and dental phenotypes.

Like most teleost fishes, cichlids possess two toothed jaws (Schaeffer and Rosen 1961; fig. 1). Cichlids have oral jaws that are homologous to mammalian jaws and that are primarily used to capture prey. Additionally, they possess toothed pharyngeal jaws that are modified gill arches (GAs) that have lost gills and are used instead to process prey (Liem 1973). In cichlids, the pharyngeal jaws are thought to be an evolutionary novelty (Liem 1973) that has been modified to be exceptionally efficient at processing food items and allowed the oral jaws to be freed up to diversify independently (Streelman et al. 2003; Hulsey 2006; Hulsey et al. 2008). The independent diversification of these two jaws has likely contributed to the outstanding diversity in shape, size, and number of cichlid teeth (Fryer and Iles 1972; Hulsey et al. 2008; Burress 2015). Therefore, the repeated evolution of similar tooth phenotypes might commonly be due to expression differences in genes that are isolated to a single jaw. Alternatively, since these two jaws are formed from serially homologous pharyngeal arches, adaptations involving their teeth might largely be generated via expression differences in genes that contribute to tooth formation on both jaws (Jernvall et al. 2000; Vonk et al. 2008; Fraser et al. 2009). For instance, the evolution of larger teeth on the pharyngeal jaw might result simply from the increased expression of genes that are highly expressed in teeth found on both jaws. The serially homologous phenotypic duplication that the cichlid pharyngeal arches offer when coupled with the rampant evolutionary replication of tooth phenotypes in cichlids provides a powerful framework for teasing apart which components of tooth development contribute to convergent adaptations.

Taking advantage of the replicated nature of cichlid diversity, we used a combination of gene expression studies and comparative genomics to identify whether there is a shared developmental genetic basis to their convergence in larger tooth sizes (fig. 1A and B). We asked several questions: 1) Do both known odontogenic genes as well as putatively new tooth genes show increased expression in the toothed oral and pharyngeal jaws as compared with the nontoothed but serially homologous GAs? 2) Are there lineage-specific genes involved in odontogenesis in cichlid fishes?, and 3) Do particular genes that are highly expressed in both toothed jaws contribute to the developmental genetic basis of convergence in pharyngeal jaw tooth size across the evolutionary diversity of cichlids? To identify genes that likely have high expression levels in cichlid teeth, we performed a transcriptomic analysis that took advantage of the serially homologous nature of cichlid tooth-bearing arches (fig. 1A). Using RNA-seq, we identified genes that are upregulated on both jaws as compared with the toothless fourth GA located posterior to the oral jaws but anterior to the pharyngeal jaws. We then examined the role of several select candidate genes during odontogenesis using in situ hybridization to validate whether these candidate genes were expressed in developing teeth. Finally, using evolutionarily disparate species pairs that differ in their pharyngeal jaw tooth size (fig. 1B), we determined if there was convergence in the levels of gene expression of two conserved and two genes present only in cichlids that formed a physical cluster in the cichlid genome.

Results

Expression of Conserved and Novel Tooth Genes in Cichlid Fishes

Based on analyses of gene expression of the toothed jaws and a homologous GA, we identified 212 genes that are differentially expressed between the lower oral jaw (LOJ) and fourth GA, of which 79 were upregulated in the LOJ. In the comparison between the lower pharyngeal jaw (LPJ) and GA, 1,099 genes were differentially expressed and 309 were upregulated in the LPJ (supplementary fig. 1, Supplementary Material online). In the LOJ, genes with functions in calcium binding and structural activity were overrepresented, but not significantly. Genes that were upregulated in LPJ compared with GA were involved in ossification, anatomical structure development, and collagen organization (supplementary fig. 2, Supplementary Material online). In total, 105 genes consistently showed lower expression, whereas 27 genes (supplementary table 1, Supplementary Material online) consistently showed higher expression in both the toothed LOJ and LPJ tissues compared with nontoothed GA (supplementary fig. 1B, Supplementary Material online). From this set of 27 genes, we identified five that were previously described to be involved in tooth development, 15 that were described to have documented functions unrelated to teeth, and seven genes that had no previous known function (supplementary table 1, Supplementary Material online).

Identification of Cichlid-Specific Tooth Genes

From the 27 genes that were upregulated in both toothed jaws, we selected four genes to further investigate based on their high expression levels in both jaws, their high fold change compared with the nontoothed GA, and their genomic location. We found that these four genes formed a compact physically linked cluster confined to a region spanning only ∼290 kb of the Oreochromis niloticus genome (fig. 3A). We examined: 1) Resuscitation promoting factor A (rpfA), 2) odontogenic ameloblast-associated protein precursor (odam), 3) secretory calcium-binding phosphoprotein 5 (scpp5), and 4) one undescribed gene (unk; supplementary table 1, Supplementary Material online). Odam and scpp5 are known to be involved in tooth formation (Kawasaki and Weiss 2003; Kawasaki et al. 2005; Lee et al. 2010, 2012). Neither unk nor rpfA, although annotated as protein coding in cichlid genomes, have been characterized previously to be involved in tooth formation or any other biological process.

We used BLAST to compare genomic and amino acid sequences of unk and rpfA from O. niloticus to other available teleost genomes and did not detect any orthologous sequences outside of cichlids. However, all cichlid species with publicly available genomes (ensembl 97, Zerbino et al. 2018) did contain orthologs of unk as well as rpfA in the same orthologous cluster with scpp5 and odam. In figure 3A, general patterns of divergence in synteny of the four genes is depicted for two evolutionarily disparate cichlid species that have well-annotated genomes, Amphilophus citrinellus from Nicaragua and O. niloticus from Africa, as well as for the damselfish Acanthochromis polyacanthus (ASM210954v1, NCBI accession number: GCA_002109545.1) and the Japanese medaka Oryzias latipes (ASM223471v1, NCBI accession number: GCA_002234715.1). The genome of the damselfish does not contain orthologous coding regions to either unk or rpfA and based on reciprocal BLASTs it possesses only small noncoding sequences (<50 bp) in this or any other region that produce hits to the two novel tooth genes (fig. 3A). Similarly, in the Japanese medaka, we could only find short noncoding orthologous sequences (fig. 3A). Additionally, in many other teleosts, this putative cichlid tooth cluster appears to show rearrangements as the region spanning unk and odam is inverted in the damselfish and many other species when compared with the genomes found in medaka and cichlids (fig. 3A). Extensive genomic comparisons therefore suggest that this cluster contains two cichlid-specific genes in a region exhibiting structural rearrangements compared with many other teleosts (see supplementary fig. 6, Supplementary Material online, for details on all investigated species).

Confirmation That Candidates Are Involved in Tooth Formation

Our in situ hybridization results confirmed expression of our candidate genes in developing cichlid teeth on both the LOJ and the LPJ (fig. 2C). The previously described tooth genes scpp5 and odam were expressed in the pulp and dental epithelium of forming teeth on both jaws. A similar pattern was observed for unk and rpfA, where both were found to be expressed in dental epithelium either lateral to the forming tooth (rpfA) or at its tip (unk; fig. 2C). Expression of these four genes was not detected in any other part of either jaw, supporting our inference that these genes are involved in tooth formation (fig. 2C).

Convergent Molecular Basis of Tooth Size

Using quantitative real-time PCR (qPCR), we analyzed expression levels of the four clustered genes (unk, odam, rpfA, scpp5) in the LPJs of 12 species that were divergent in pharyngeal tooth size. In each of the six species pairs, one species with enlarged molar-like teeth and one species with smaller pointed teeth were examined (fig. 3B), and their pharyngeal tooth sizes were also confirmed to be significantly different (fig. 3B). Expression levels of the four genes commonly showed significant differences (P < 0.01) among members of the species pairs (supplementary table 4, Supplementary Material online). For unk, expression differences occurred in the opposite directions for some species pairs. But, we found a pattern for odam in which large-toothed species either showed statistically indistinguishable or higher expression compared with smaller toothed species. However, across all six species pairs, rpfA and scpp5 showed higher expression in the larger toothed species (fig. 3B). This gene cluster contributes to a convergent developmental genetic basis for increased tooth size across cichlid diversity.

Discussion

Developmental diversity of teeth in cichlids and likely other nonmammalian vertebrates is probably greater than previously appreciated. Our comparisons of the genes in the serially homologous pharyngeal elements of cichlids allowed us to identify a substantial number of known odontogenic and putatively new tooth genes that could be involved in forming cichlid teeth. Despite being expressed in both toothed jaws, some of these novel genes appear to be exclusively found in cichlid genomes and were confirmed to be deployed extensively during odontogenesis. A genomic cluster containing a combination of two of these novel genes as well as highly conserved odontogenic genes were found to have gene expression changes that are associated with convergent evolution of cichlid tooth-size diversity. The evolutionary coupling of the cichlid pharyngeal jaw’s anatomical novelties with this “tooth size” genomic cluster could have provided a key mechanism for cichlid’s remarkable evolutionary success.

The genes we identified as being highly expressed in toothed jaws compared with nontoothed GAs provide evidence that many previously unidentified genes could be significantly involved in cichlid tooth formation (supplementary table 1, Supplementary Material online). Although they are embedded within much less dynamic and mineralized jaw bones, formation of replacement teeth in cichlids involves production of various tissues that are highly transcriptionally active (Rasch et al. 2016). However, not all of these genes could be involved in odontogenesis as there are structures other than teeth that these highly expressed genes could contribute to in the jaws (Alappat et al. 2003). Nevertheless, our stringent filtering probably led to an underestimation of the number of genes involved in odontogenesis. If these expression results are as generalizable to other toothed vertebrate structures as we suspect, then our transcriptomic analyses suggest that evolutionarily and phenotypically distinct mammalian models of tooth development could represent only a small fraction of the developmental genetic diversity of vertebrate teeth. This finding is likely to be true for many traits such as paired limbs (Stopper and Wagner 2005), chambered hearts (Olson 2006), and teeth (Kawasaki et al. 2004) that despite having diversified extensively likely had a single origin during vertebrate evolution. In general, the diversity of genes expressed during organogenesis for many characteristically vertebrate structures might be far greater than a biased focus on mostly mammalian models of development has led us to believe.

Our results provide a putative example where lineage-specific developmental mechanisms were integrated into otherwise highly conserved developmental programs. Based on their physical linkage and shared high expression levels in transcriptomes from both cichlid jaws, unk and rpfA might be inferred to be coregulated with the more conserved and more widely deployed tooth genes, scpp5 and odam (Kawasaki et al. 2005; Lee et al. 2010). Further, one might infer from the function and location of expression of numerous other tooth genes (Kawasaki et al. 2005) that unk and rpfA could be similarly involved in mineralization of enameloid or production of the enameloid organic base (Kawasaki et al. 2005; Assaraf-Weill et al. 2014). However, protein-coding regions of both genes lack calcium-binding motifs characteristic of genes involved in tooth mineralization (Kawasaki and Weiss 2003). In addition, our in situ hybridization results do show that the novel tooth genes exhibit somewhat unique patterns of expression in developing teeth that are distinct from scpp5 and odam (fig. 2C). Additionally, their minor comparative discrepancies in relative expression across species pairs in comparisons with the two more conserved genes suggests the two novel genes likely are subjected to some independent regulation (fig. 3B). The relative independence of gene regulation is known to be critical to whether gene clusters promote adaptation (Singer et al. 2005; Lotterhos et al. 2018; Zeng et al. 2018). However, our results provide evidence that genomic clustering as observed for these four physically linked genes could be a common mechanism through which novel genes are incorporated into conserved but adaptive developmental programs like those for tooth morphogenesis.

How novel genes like unk and rpfA originate is often unclear, and we can envision multiple possible scenarios. Gene duplication with subsequent neofunctionalization is a frequently observed means of creating novel genes (Ohno 1970). But, we did not identify any paralogous sequences of unk or rpfA in cichlid genomes or any orthologous sequences in published genomes outside of the cichlid lineage (fig. 3A). It is also possible that unk and rpfA could have originated from noncoding genomic elements in the ancestor of all cichlids (Wu et al. 2011; Neme and Tautz 2013), but the existence of several exons in each gene would seem to argue against this. As the genomes of additional teleosts that are more closely related to cichlids become available, we might be able to clarify the origin of these genes. Additionally, these are not likely the only novel genes that contribute to cichlid diversification. The birth of new genes was once regarded as a rare event and generally was inferred to have little influence on how phenotypes evolve (Jacob 1977). But, in recent years evidence has accumulated that gene birth occurs more frequently than previously thought (Wu and Zhang 2013) and these recently originated genes can rapidly acquire functions that contribute to adaptation (Chen et al. 2010; Wu et al. 2011; McLysaght and Guerzoni 2015). The unk and rpfA genes appear to commonly contribute to the repeated enlargement of teeth that permit the processing of more durable prey items. Thus, these types of novel lineage-specific genes might generally play an important role in organismal diversification and in convergently evolved adaptations.

The pharyngeal jaw that bears the repeatedly evolved enlarged teeth examined here is often considered to be a key innovation that facilitated the immense trophic diversity and even speciation of cichlids (Liem 1973; Kocher 2004; Hulsey et al. 2008). Expression of genes in the cluster is associated with tooth size divergence and could have evolutionarily augmented the macroscopic anatomical specializations that constitute this key innovation. The novel insertion of the pharyngeal muscles, unique upper jaw joints, and sutured medial fusion of the LPJ all contribute to the ability of cichlids to crush durable prey, an ability that is rare in teleosts (Liem 1973). The “tooth size” genomic cluster investigated here therefore might have repeatedly and synergistically interacted with these other characteristic cichlid traits and facilitated the evolutionary likelihood that these fishes will readily re-evolve enlarged teeth that are commonly associated with a crushing trophic habit. The cichlid pharyngeal jaw has been proposed to enhance the evolvability, or capacity for adaptive evolution, of this family of teleosts and this gene cluster containing cichlid-specific genes could have amplified this capacity (Galis and Metz 1998). As we continue to genomically dissect the developmental mechanisms generating convergent adaptations that characterize the exceptional diversity and species richness of lineages such as cichlid fishes, it will be key to more extensively recognize the genomic and developmental interplay between evolutionary conservation and novelty.

Materials and Methods

All fish were raised in aquaria at the University of Konstanz animal care facility. Prior to tissue extraction, specimens were euthanized with an overdose of tricain (MS-222). For RNA-sequencing and in situ hybridization, we used the African cichlids Astatoreochromis alluaudi and Astatotilapia burtoni, respectively, and species listed in figure 1B were used for quantitative PCR analysis.

Gene Expression Analyses

For RNA-seq analysis, we dissected three distinct tissues from Astatoreochromis alluaudi including the LOJs (n = 10), the most posterior GA (n = 10) and LPJs (n = 20) and stored them in RNAlater (Sigma-Aldrich, St. Louis). These fish ranged in size from 46 to 110 mm. Because cichlids replace their teeth throughout their lifetimes, teeth of fish at these sizes should show high transcriptional activity and allow us to primarily observe expression in tooth genes. Prior to storage of tissue in RNA-later at −20 °C, we removed soft connective and muscle tissue from the elements. We then extracted RNA using the ReliaPrep miRNA cell and tissue miniprep system (Promega, Madison, WI) according to Schneider et al. (2014). We then prepared RNA libraries for each sample with the SENSE mRNA library preparation kit (Lexogen, Vienna, Austria) with an input of 100 ng total RNA per sample and unique barcode sequences for each individual sample. After equimolar pooling of RNA libraries, we size-selected all fragments in a range of 300–600 bp and sequenced libraries (paired-end 2× 150 bp) on an Illumina Hiseq X-ten platform at the Beijing Genomics Institute (BGI, Hong Kong, China).

For quality control, we removed low-quality bases from each read. Bases at the beginning and end of each read with a quality score <3 were removed. Each read then was scanned with a four-base sliding window and bases were clipped if the average quality score within the sliding window dropped <15. Further, adapter sequences were clipped from each read. Afterward, all reads shorter than 145 bp were removed. Remaining high-quality reads were mapped to the publicly available O. niloticus reference transcriptome (Orenil_1.0, GCA_000188253.1) using the program Salmon v0.8.1 (Patro et al. 2017) that allowed us to obtain normalized expression values for each transcript. For differential expression analyses, we used the Bioconductor package limma (Ritchie et al. 2015) in the R statistical environment and discarded loci falling under a minimum threshold of 20 TPM in at least five samples per tissue. The high threshold was used to focus on highly expressed genes and to reduce the analyzed set of genes from 33,220 to 6,407. Differential expression for the remaining genes was analyzed for two contrasts: 1) genes differentially expressed between the LOJ and GA as well as 2) genes differentially expressed between LPJ and GA (fig. 1A). We concentrated on these contrasts involving the GA because this structure represents a nontoothed serially homologous element that shares developmental similarity with both of the toothed jaws. Gene-wise linear models with tissues as explanatory variables and expression of each gene as the response variable, as implemented in limma, were used to determine significance over all three tissues to allow for integrated inferences over the whole experiment (Law et al. 2014). P values were further corrected using the Benjamini–Hochberg method (Benjamini and Hochberg 1995) (FDR corrected P values below α = 0.05). An absolute log₂-fold change threshold (|log₂ = 1.5|) was applied on significantly differentially expressed genes for each independent contrast to focus on genes highly differentiated among the toothed and nontoothed structures. Gene ontology enrichment analysis was performed as described in Hart et al. (2018) using g: Profiler (version e98_eg45_p14_bca6d38, Raudvere et al. 2019). We used the g: SCS method, taking into account the structure of functional annotations, to determine significant terms with P value <0.05 (Reimand et al. 2007). Gene functions were visualized with REVIGO (Supek et al. 2011).

Comparative Genomics

The comparative analysis of the investigated genomic cluster containing the four genes unk, odam, rpfA, and scpp5 was performed using BLAST (2.9.0) for similar (megablast) as well as more dissimilar sequences (discontinuous megablast). We compared the genomes of 23 different fish species including six cichlid species (supplementary fig. 6, Supplementary Material online) against which we blasted the whole genomic region, coding sequences, individual exons as well as amino acid sequences from O. niloticus, respectively. We only considered hits significant if they were >50 bp with an E-value <1.0⁻¹⁰. Genomic locations for each gene were extracted from publicly available annotations in the ensemble database (Zerbino et al. 2018; ensembl 97).

Fluorescent In Situ Hybridization

We performed in situ hybridization of putative tooth genes on cryo-sectioned LOJ and LPJ of adult Astatotilapia burtoni (n = 6). Jaws were fixed in 4% paraformaldehyde in phosphate-buffered saline (PBS) overnight at +4 °C. Following fixation, jaws were demineralized in 10% acetic acid in 4% PFA in PBS for 7 days at 4 °C. Jaws were then incubated in 18% sucrose in PBS overnight at 4 °C to increase tissue stability during freezing and sectioning. Afterward, we incubated jaws in 50% Tissue-Tek O.C.T. Compound (Sakura Finetek, The Netherlands) for several hours prior to embedding the jaws in pure Tissue-Tek O.C.T. Compound at −80 °C. We prepared 25 µm cryo-sections of the jaw tissues at −20 °C and stored the sections until further use at −80 °C.

Jaw sections were stained for either one of the four closely examined genes (unk, odam, rpfA, scpp5) with a fluorescently labeled probe according to the protocol in Woltering et al. (2009). Probes were labeled with DIG-labeled UTP during reverse transcription of template DNA. Hybridization was performed at 65 °C overnight, and subsequent antibody incubation was performed overnight at room temperature with an Anti-Digoxigenin-Peroxidase (Merck, Germany) for detection of digoxigenin-labeled probes. Fluorescent detection was performed with a thyramid amplification system (TSA Amplification System, Perkin Elmer) for 20 min at room temperature. Both bright-field and fluorescent photographs of in situ gene expression were obtained with a Leica DM6B microscope and original colors replaced in ImageJ (Schindelin et al. 2012) for presentation in figure 2.

qPCR Assay for Validation of Differential Gene Expression

To determine if large-toothed pharyngeal jaws convergently differed in gene expression levels, we examined 12 different cichlid species that form six pairs of closely related species. These pairs were chosen to be phylogenetically distributed across cichlid diversity (Friedman et al. 2013). The pairs of Amphilophus citrinellus and Archocentrus centrarchus, Amatitlania siquia and Hypsophrys nematopus, and Thorichthys meeki and Trichromis salvini are native to the Neotropics. Neolamprologus christy and N. brichardi, Astatoreochromis calliptera and Melanochromis auratus, and Mylochromis mola and Copadichromis borleyi are all native to Africa. Six individuals per species were examined except for the species Thorichthys meeki, Trichromis salvini, N. brichardi, and M. mola in which five individuals were examined. Pairs were chosen specifically because previous observations indicated that one species likely exhibited enlarged molar-like teeth and the other displayed only smaller, more pencil-like teeth (see full list of species in fig. 3B) on the pharyngeal jaw. All individuals used were adult (>50 mm) lab-reared fish. We photographed each fish as well as their respective LPJ from a dorsal view using a Leica MZ10F stereomicroscope. These digital photographs were imported into Image J (Schindelin et al. 2012) and the standard length as well as dorsal tooth areas was measured. The area of the six most posterior teeth along the midline of the jaw that tend to be the largest teeth and are likely to be homologous across all species were averaged for each individual. These average tooth areas were then square root transformed and adjusted for size by taking residuals of the regression on individual standard length for each species pair. Subsequently, we extracted RNA from LPJs as described above (see Gene Expression Analyses).

Extracted RNA was reverse transcribed using the GoScript Reverse transcription system (Promega, Madison, WI) with an input of 200 ng total RNA according to manufacturer’s protocol. Synthesized cDNA was diluted 1:5 in nuclease free water and 2 µl were used for qPCRs according to manufacturer’s protocol (GoTaq qPCR MasterMix, Promega, Madison, WI). Expression levels of the four candidate genes (unk, odam, rpfA, scpp5) as well as two housekeeping genes (actinR, twinfilin; Gunter et al. 2013) were analyzed by qPCR with the GoTaq qPCR System (Promega, Madison, WI) as described in Karagic et al. (2018). Significant differences in each genes’ expression values were compared for each species pair using individuals as replicates with a t-test (α = 0.05). Since analysis of replicate sister taxa provides one of the most robust frameworks for phylogenetic comparative methods (Felsenstein 1985), we also generated cichlid-wide P values over all species pairs globally with an ANOVA to assess if gene expression showed cichlid-wide evidence of convergent associations with tooth size.

Supplementary Material

Supplementary data are available at Molecular Biology and Evolution online.

Acknowledgments

We thank the International Max-Planck Research School for Organismal Biology, the helpful staff of the animal care facility at the University of Konstanz and Paolo Franchini for valuable comments. This work was supported by the Deutsche Forschungsgemeinschaft (DFG HU2278/1-1 to C.D.H. and DFG ME1725/20-1 to A.M.) and the European Research Council (ERC advanced grant No. 293700 to A.M.).

References

Author notes