Skeletal Mineralization in Association with Type X Collagen Expression Is an Ancestral Feature for Jawed Vertebrates

Abstract In order to characterize the molecular bases of mineralizing cell evolution, we targeted type X collagen, a nonfibrillar network forming collagen encoded by the Col10a1 gene. It is involved in the process of endochondral ossification in ray-finned fishes and tetrapods (Osteichthyes), but until now unknown in cartilaginous fishes (Chondrichthyes). We show that holocephalans and elasmobranchs have respectively five and six tandemly duplicated Col10a1 gene copies that display conserved genomic synteny with osteichthyan Col10a1 genes. All Col10a1 genes in the catshark Scyliorhinus canicula are expressed in ameloblasts and/or odontoblasts of teeth and scales, during the stages of extracellular matrix protein secretion and mineralization. Only one duplicate is expressed in the endoskeletal (vertebral) mineralizing tissues. We also show that the expression of type X collagen is present in teeth of two osteichthyans, the zebrafish Danio rerio and the western clawed frog Xenopus tropicalis, indicating an ancestral jawed vertebrate involvement of type X collagen in odontode formation. Our findings push the origin of Col10a1 gene prior to the divergence of osteichthyans and chondrichthyans, and demonstrate its ancestral association with mineralization of both the odontode skeleton and the endoskeleton.


Introduction
Several factors make type X collagen an excellent candidate to study the evolution of the genetic mechanisms involved in skeletogenesis. First, although type X collagen was first identified more than 30 years ago in chicken long bones (Schmid and Linsenmayer 1985), it has to date only been described in osteichthyans. An immunohistochemical analysis using rabbit anti-Col X antibodies labeled the shark skeleton (Seidel et al. 2017), but the elephant shark and whale shark genome sequencing projects have failed to identify type X collagen genes (Venkatesh et al. 2014;Read et al. 2017). Second, type X collagen is a crucial factor of skeletal mineralization in amniotes, expressed in hypertrophic chondrocytes during endochondral ossification (Schmid and Linsenmayer 1985;Jacenko et al. 1993;Mcintosh et al. 1994;Kwan et al. 1997), and contributes to the initiation of calcium precipitation in the cartilaginous extracellular matrix (Kirsch and Wuthier 1994). Third, studies in nonamniote osteichthyans have revealed that Col10a1 is not limited to cartilage mineralization as it is expressed by osteoblasts in the amphibian Xenopus tropicalis (Aldea et al. 2013) and the actinopterygians: the medaka Oryzias latipes (Renn and Winkler 2010), the spotted gar Lepisosteus oculatus, and the zebrafish Danio rerio (Laue et al. 2008;Albertson et al. 2010;Eames et al. 2012). This suggests that in osteichthyans, type X collagen was originally involved in the mineralization of both the cartilage and bone matrix of the endoskeleton. Finally, Col10a1 originated through the two rounds of genome duplication that occurred early in vertebrate evolution, together with its paralogous Col8a1 and Col8a2 genes (Aldea et al. 2013). This suggests that Col10a1 should have been ancestrally present in chondrichthyans, and raises the question of when it was recruited to the skeletal mineralization processes.
Chondrichthyans or cartilaginous fishes (Elasmobranchii and Holocephalii) have lost perichondral bone and all large dermal bones (Janvier 1996). Chondrichthyan skeletons, however, are characterized by several types of mineralized tissues: globular mineralized cartilage also found in fossil jawless vertebrates (Ørvig 1951;Dean and Summers 2006) and prismatic cartilage which is a chondrichthyan-derived feature (Seidel et al. 2016). A third type of mineralization is only found in elasmobranchs (sharks, rays, and skates), in the fibrous tissue surrounding their notochord, the site of areolar mineralization (Arratia et al. 2001;Eames et al. 2007;Enault et al. 2015;Criswell et al. 2017;Atake et al. 2019). Finally, a fourth mineralized tissue was described in the neural arches of few elasmobranch species: it resembles lamellar bone, and will be here-in named "lamellar mineralization" as its homology with bone in osteichthyans is still debated (Wurmbach 1932;Lallier et al. 1982;Peignoux-Deville et al. 1989;Eames et al. 2007;Enault et al. 2015). Nevertheless, the site of densest mineralization in past and extant vertebrates is found in their odontode skeleton, which includes teeth, scales, and fin spines (Qu et al. 2015). These structures are made of enamel or enameloid outer layer and a dentin core produced as a result of epithelial-mesenchymal interaction, where ameloblasts of epithelial layer secrete enamel matrix and odontoblasts of mesenchymal layer secrete dentin matrix (Hall 2015).
Here, we report the identification of several type X collagen genes in both clades of cartilaginous fishes, elasmobranchs (including the lesser spotted catshark Scyliorhinus canicula, whale shark Rhincodon typus, thornback ray Raja clavata, little skate Leucoraja erinacea, thorny skate Amblyraja radiata) and holocephalans (represented by the elephant shark Callorhinchus milii). We conclude that these type X collagen genes arose from a series of local tandem duplications that occurred prior to the divergence of the holocephalan and elasmobranch lineages, in the time period when the unique hard-tissue repertoire of chondrichthyans evolved. The major expression site of these type X collagen genes in the cartilaginous fishes is in the teeth and scales of the odontode skeleton, whereas only one gene duplicate is also expressed in vertebral mineralizing tissues of the endoskeleton. We further show that a type X collagen expression site is retained in teeth of two osteichthyans, the zebrafish D. rerio and the western clawed frog X. tropicalis. Our results indicate an ancestral shared involvement of Col10a1 genes in the formation of a mineralized odontode skeleton, in addition to its role in endoskeletal mineralization, as early as the last common ancestor of all jawed vertebrates.

Duplication of Col10a1 Gene in Holocephalan and Elasmobranchs
We searched the whole elephant shark genome assembly and identified five Col10a1 genes on scaffold 71, one Col8a1 gene on scaffold 42, and one Col8a2 gene on scaffold 121. Using these predicted coding sequences (CDS) to screen transcriptomic data of other chondrichthyan species led to the identification of six similar coding sequences in the catshark S. canicula and thornback ray Ra. clavata, as well as a single Col8a1 and Col8a2 sequence. In the transcriptome of little skate L. erinacea we have identified partial gene transcript sequences LSb2-ctg50203 for Col10a1.1; LSb2-ctg2668069 for Col10a1.2; LSb2-ctg41083 for Col10a1.4; LSb2-ctg68055 for Col10a1.5; LSb2-ctg688724 for Col10a1.6; LSb2-ctg44808 for Col8a1; and LSb2-ctg1616013 for Col8a2.
In whale shark genome assembly we identified one of the hits, which was a long predicted transcript (XM_020514194.1) spreading over 100 kb on the genomic locus NW_018064585.1, to which four of the catshark cDNAs could be aligned in tandem (fig. 1 and supplementary table 1, Supplementary Material  online). Two additional predicted CDS had a positive reciprocal BLAST best hit with catshark Col10a1 copies on two other genomic loci ( fig. 1). Together, we identified six Col10a1, one Col8a1, and one Col8a2 gene sequences in the whale shark genome.

Conserved Synteny of Col10a1 Genes in Jawed Vertebrates
The human genomic locus of the Col10a1 gene (Ensembl ENSG00000123500) is located in the 103 kb long intron 6 of the 5 0 -nucleotidase domain containing 1 gene (Nt5dc1, Ensembl ENSG00000178425). These two genes are located on the opposite strands and a Fyn related Src family tyrosine kinase gene (Frk ENSG00000111816) has the same orientation as Col10a1 but located 40 kb away from the 5 0 end of the Nt5dc1 gene ( fig. 1). Conserved synteny is described in the Xenopus, gar, and coelacanth genomes ( fig. 1). In the zebrafish D. rerio, the third round of teleost-specific genome duplication led to two Col10a1 copies located on chromosomes 17 (col10a1a) and 20 (col10a1b), the first one located 14 kb upstream of one copy of Frk (fyna ENSDARG00000011370), the latter is located in the intron 6 of nt5dc1 (ENSDARG00000006797) and 20 kb upstream of another copy of Frk (fynrk ENSDARG00000027807).
We first focused on intron 6 of Nt5dc1 in the elephant shark C. milii to undergo the Genscan analysis that led to the identification of five consecutive Col10a1 genes ( fig. 1). We further analyzed the 85 kb long intergenic region between the first exon of Frk and the first exon of Nt5dc1 but found no Col10a1 copy (this region however also contains a sequencing gap). In the whale shark R. typus, strong similarity was detected in several predicted transcripts: LOC109931113; LOC109919892 located 30 kb upstream of a Frk gene (Gene ID: 109919891); LOC109914588 which is a composite predicted transcript including four consecutive Col10a1 sequences, and is located in the intron 6 of Nt5dc1 (LOC109914587, fig. 1). The cDNA sequences from the catshark transcriptome were used as BLAST queries against a partial catshark genomic assembly (Oulion et al. 2010) and could be located on different short scaffolds. One scaffold confirmed the synteny of one copy with the six first Nt5dc1 exons, and another scaffold showed the tandem position of the Col10a1.3 and Col10a1.4 copies ( fig. 1). As no fully assembled batoid genome is currently available, and to better support synteny of the tandem copies in elasmobranchs, we used cDNA sequences from the thornback ray and little skate as BLAST queries against the draft genome sequence of the thorny skate Amblyraja radiata (assembly version sAmbRad1_p1): exons of all six tandem gene copies of Col10a1 could be identified along a 400 kb long portion of the scaffold 000143F in synteny with Frk and Nt5dc1 genes ( fig. 1). The data obtained in whale shark and thorny skate genomes highlights the position of the Col10a1.1 copy, not in the intron 6 of Nt5dc1 but in the intergenic space between the start of genes Frk and Nt5dc1 ( fig. 1). Genomic data for osteichthyan species were extracted from public databases (Ensembl, NCBI, XenBase) and analyzed to reconstruct the synteny blocks of duplicated Col10a1 genes in chondrichthyans. Distances are to scale within each species but may differ between species. Abbreviation: ex, exon.
Skeletal Mineralization in Association with Type X Collagen . doi:10.1093/molbev/msz145 MBE Col10a1 group was recovered monophyletic with strong support. Within it, all osteichthyan Col10a1 sequences grouped together and were sister-group to the monophyletic group of the chondrichthyan Col10a1 sequences. The branching between several paralogous groups suggests at least five local duplication events that most probably occurred in the chondrichthyan stem lineage as four out of six clades include both elasmobranch and holocephalan sequences. We chose to name the tandem gene duplicates by adding a number, from .1 to .8, following the original gene name Col10a1, with respect to the time of divergence (following the phylogenetic reconstruction [ fig. 2] and their respective position along the chromosome [ fig. 1]). The first duplication led to the Col10a1.1 copy located in the Nt5dc1-Frk intergenic sequence and later lost in holocephalans, whereas its sister copy remained in the intron 6 of Nt5dc1

Protein Domain Organization of Col10a1 Duplicates in Chondrichthyans
Protein domain predictions showed that the elephant shark and catshark Col10a1 duplicates share typical domain architecture and organization with osteichthyan Col10a1 proteins ( fig. 2). Missing signal peptides in some sequences can be explained by incompleteness of genomic and cDNA data. In both holocephalans and elasmobranchs we detected an additional domain in Col10a1.5, located after the signal peptide and before the triple alpha helix repeat ( fig. 2). This domain (named SE) is rich in serine (S), glutamic acid (E), and also to some degree in aspartic acid (D) (supplementary fig. 3, Supplementary Material online). Conserved domain analysis annotated this domain as DMP1 superfamily cl25845 and pfam07263 (Marchler-Bauer et al. 2017). This superfamily includes dentin matrix acidic phosphoprotein 1 (DMP1) which has a similar composition of amino acids. For example human DMP1 has the highest number of S (21%), E (15%), and D (13%). In comparison, SE domains in Col10a1.5 of catshark, whale shark, thornback ray, and elephant shark contain the highest number of S (21-31%) and E (10-12%) but lower number of D (5-8%) (

Specialized Expression Patterns of Duplicated Col10a1 Genes in Elasmobranchs
From the transcriptome data of the catshark lower jaw, the catshark vertebra and thornback ray lower jaw, we presented Col10a1 gene expression values as Transcripts Per Million (TPM). In lower jaws all six tandem duplicates were expressed, with values above 10 detected for Col10a1.1, Col10a1.2, Col10a1.4, and Col10a1.6 in catshark and ray, as well as for Col10a1.3 in catshark (table 1). In contrast, in catshark vertebra only Col10a1.4 had a value above 10, although Col10a1.2 also showed detectable levels below 10.
To further locate the detailed gene expression pattern of each duplicate, we performed in situ hybridization with RNA probes specific to each copy of catshark Col10a1 MBE duplicate could not be detected in vertebral tissues (data not shown) despite low TPM values for this gene in the vertebral column transcriptome. We consider that little amounts of skin tissue including scales in RNA-seq samples might be the reason for these TPM results. The sites of mineralization described in catshark vertebrae (Enault et al. 2015) include the ring of fibrocyte-like cells, involved in areolar mineralization surrounding the notochord, the fibrous tissues surrounding the cartilaginous neural arches and the cartilaginous core of the neural arches, not mineralized at that stage but known to Bipartition support was evaluated by 100 bootstrap replicates. Maximum support is shown by filled circles. Bootstrap supports below 75 are not displayed, and the branch thickness of their corresponding bipartition was decreased. Col8a1 and Col8a2 genes were used as outgroups and tree was rooted following Aldea et al. (2013). Predicted protein domains for catshark, elephant shark, and human Col10a1 sequences are displayed next to the gene names. Protein domains are labeled as SP-signal peptide, GXY-collagen triple alpha helix repeat, C1q-complement component 1q domain, and SE-serine-rich, acidic domain. Abbreviations: C., Chondrichthyes (cartilaginous fishes); O., Osteichthyes (bony fishes).
Skeletal Mineralization in Association with Type X Collagen . doi:10.1093/molbev/msz145 MBE and col10a1b during larvae development. Expression of col10a1a was detected in bones (e.g., dentary, maxilla, ethmoid plate, opercle, and cleithrum) of the zebrafish larvae ( fig. 4A and B). Additional expression sites were detected in pharyngeal teeth at days 4, 5, and 6 of zebrafish development figure 4A-C. Expression of col10a1b was detected only in teeth at days 4 and 6 ( fig. 4D-F). Post in situ paraffin sections reveal the expression of col10a1a and col10a1b in the mesenchyme of the older tooth (with deposited matrix) and the expression of col10a1b in the epithelium of a less developed tooth bud (without deposited matrix, fig. 4C and F).

Discussion
Our analysis of genome and transcriptome data from a variety of chondrichthyan species led to the identification of five Col10a1 gene duplicates in holocephalans (elephant shark), and six Col10a1 gene duplicates in elasmobranchs (lesser spotted catshark, whale shark, thornback ray, little skate, and thorny skate). This is a first reported characterization of the Col10a1 gene orthologues outside of the osteichthyan group, which supports the origin of Col10a1 before the divergence of osteichthyans and chondrichthyans (dated to 422-463 Ma, Benton et al. 2009). In addition, synteny and phylogenetic data (figs. 1 and 2) support an early timing of these duplications, before the divergence between holocephalans and elasmobranchs (dated to 352-447 Ma, Benton et al. 2009), with little subsequent rearrangement.
Our data suggest that Col10a1.1 and Col10a1.4 might have later been lost in the holocephalan lineage, although they may be missing simply because of incomplete genome sequencing. A putative secondary duplication of Col10a1.3 might also be interpreted from our synteny data and phylogenetic reconstruction, with elephant shark Col10a1.7 and Col10a1.8 duplicates clustering in a sister-group to the elasmobranch Col10a1.3 sequences. However this topology is poorly supported with nucleotide sequence, and not recovered in the analysis of the protein sequences (supplementary fig. 1, Supplementary Material online) so we chose to name them MBE independently. The identified chondrichthyan duplications mostly occurred as tandem duplications at the locus conserved for the Col10a1 gene in bony fishes, except for the first event leading to the Col10a1.1 copy, which is found outside on the Ntd5dc1 intron 6, however still located in synteny, close to the Frk gene ( fig. 1 and supplementary fig. 2, Supplementary Material online). The initial conservation of such tandem duplicates might be explained by a positive selection for a higher copy number (Hahn 2009;Schrider and Hahn 2010;Margres et al. 2017) given the importance of high protein quantity in the synthesis of extracellular matrix to be provided during development. However, based on our expression data in shark, the long-term conservation of these copies seems to be due to subsequent subfunctionalization including the acquisition of specialized expression patterns and/or protein function. The expression of all six catshark Col10a1 gene duplicates in the developing teeth and scales first brings up their importance for the chondrichthyan odontode skeleton development, greater than their involvement in endoskeletal genesis. We discovered several levels of potential functional diversification represented by epithelial, epithelial-mesenchymal, and mesenchymal gene expression. Three duplicates Col10a1.1, .2, and .3 showed transcription in ameloblasts of the catshark teeth and scales: these cells are secretory and/or maturation (mineralizing) ameloblasts (Kemp 1985). Previous evidence supports the ectodermal origin of the mineralized components of the shark tooth enameloid however there is an ongoing discussion whether organic component is secreted by ameloblasts or odontoblasts, or both (Kemp 1985). Taken together, our data suggest the involvement of catshark Col10a1.1, Col10a1.2, and Col10a1.3 in the enameloid formation that had never been shown before for nonfibrillar collagens of type X.
On the other hand we describe four of the duplicated genes (Col10a1.3 to Col10a1.6) with expression in mesenchymal cells of developing teeth and scales in the catshark. Interestingly, patterns in the mesenchymal compartment showed variation-pulp of the tooth and scale for Col10a1.6 versus extensive expression in the root region for Col10a1.3 and Col10a1.4. Nevertheless, these are difficult to interpret, as the dynamics of odontoblasts and adjacent cell type differentiation have not been studied in detail in chondrichthyan fishes. Still, the expression of Col10a1.5 seems to be concentrated to the secretory odontoblasts indicating the involvement in the production of dentin matrix.
We questioned protein function for these various duplicates, by predicting protein domains, and identified collagenous domains with length variation between duplicates. In model organisms, type X collagen was linked to biomineralization processes through its ability to bind with annexins from the matrix vesicles that mediate apatite precipitation in the extracellular matrix (Kirsch and Wuthier 1994 and reviewed by Bottini et al. 2018). However the exact functions of the type X collagen protein domains in this process remain unknown and we can only speculate on these functions in chondrichthyans.
Unexpected identification of a novel serine-rich acidic domain in Col10a1.5 with similarity to DMP1 could indicate the combined role for this duplicate in the deposition of both acidic and collagenous matrix in the chondrichthyan dentin. DMP1 belongs to the acidic subfamily of secretory calcium-binding phosphoproteins (SCPP) together with bone matrix proteins (Kawasaki 2009). Strikingly, SCPPs could not be identified in chondrichthyan genomes (Venkatesh et al. 2014;Enault et al. 2018). The appearance of the SE domain might therefore come from a convergent event of protein sequence evolution that allows structural similarity between dentin tissues of osteichthyan and chondrichthyan groups. Overall, our analysis showed that these genes that are transcribed during the odontode skeleton development arose through rapid tandem duplication events and underwent subsequent subfunctionalization acquiring some specific characteristics that are currently unique for chondrichthyan species.
Analysis of Col10a1 expression data only in catshark would have suggested that an ancestral type X collagen gene had been recruited to the odontode skeleton in the chondrichthyan lineage. However, a previous study reported Col10a1 expression in zebrafish teeth (Simões et al. 2006) and our detailed analysis of Col10a1 expression in zebrafish Skeletal Mineralization in Association with Type X Collagen . doi:10.1093/molbev/msz145 MBE and Xenopus also revealed expression sites in the ameloblasts and odontoblasts of osteichthyan teeth. This supports that ameloblastic and odontoblastic expression of the Col10a1 is a shared character of bony fishes (although apparently lost in mammals) and cartilaginous fishes, most probably inherited from the last common ancestor of all jawed vertebrates ( fig. 5).
Because the osteichthyan type X collagen is a marker for osteoblasts and hypertrophic cartilage in osteichthyans, we questioned the importance of type X collagen for chondrichthyans where endochondral bone is missing. Only one duplicate (Col10a1.4) was strongly expressed during endoskeletal mineralization in the catshark. Expression was detected in standard chondrocytes going through globular mineralization, in "fibrocyte-like chondrocytes" of the fibrous vertebral body and in "osteocyte-like" cells that surround the neural arches. We did not detect expression of Col10a1.4 or any other Col10a1 duplicate in the Meckel's cartilage, which is not mineralized yet, although the initiation of mineralization is set up at that developmental stage (Enault et al. 2015). We cannot completely exclude the possibility of convergent recruitment of type X collagen in several cell lineages and species lineages over vertebrate evolution, but this hypothesis involves many events: three independent recruitments in ameloblasts, odontoblasts, and mineralizing chondrocytes in early jawless vertebrates; osteoblasts in bony fishes; vertebral body chondrocytes in elasmobranchs and osteocyte-like cells in some shark species.
Instead, we propose an alternative scenario where the sharing of type X collagen gene expression by all these mineralizing cells is a hallmark of an ancestral biomineralization toolkit in early vertebrates that have been deployed in all mineralizing cells (fig. 5). Tandem duplication of type X collagen genes in chondrichthyan lineage was consequently followed by subfunctionalization of tandem duplicates. This process most probably required diversification of cisregulatory mechanisms (Baudouin-Gonzalez et al. 2017). In this scenario, the regulatory elements linked to the ancestral Col10a1 gene were activating expression in endoskeletal mineralizing cells, as well as ameloblasts and odontoblasts. These regulatory elements were then differentially inactivated between tandem duplicates so that only Col10a1.4 kept the ancestral endoskeletal regulatory elements whereas all duplicates kept and also shared ameloblast and/or odontoblast specific regulatory elements.
As a final conclusion, we speculate here that the ancestral jawed vertebrate Col10a1 gene originally belonged to a shared gene regulatory network linked to extracellular matrix mineralization. This toolkit was then initiated in a variety of cells found in the odontode skeleton, dermal skeleton, and endoskeleton of early vertebrates. Our analysis of Col10a1 subfunctionalization after tandem duplications in chondrichthyans is a great example of how such toolkit can diversify during evolution. This scenario had never been explored before because of the lack of type X collagen expression by osteoblasts and odontoblasts in amniotes, in contrast to the situation of type I collagen which is involved in matrix secretion for dentin, bone, and even amphibian enameloid (Assaraf-Weill et al. 2014). However, type I collagen has no function in cartilaginous matrix and is never found associated to mineralized cartilage (Eames et al. 2007;Enault et al. 2015;G omez-Picos and Eames 2015;Cervantes-Diaz et al. 2017). For these reasons, candidates for the description of the ancestral gene regulatory network linked to biomineralization should be chosen for their interaction with type X collagen, and be tested in rayfinned fishes and cartilaginous fishes rather than amniotes.

Embryo Collection and Ethics Statement
Embryos of the lesser spotted catshark S. canicula from a Mediterranean population were obtained at the Station M editerran eenne de l'Environnement Littoral, Sète, France. Collection and handling of catshark embryos were carried out in full compliance of institutional, national, and international guidelines (European Communities Council Directive of 22 September 2010 [2010/63/UE]) and did not require further approval by an ethics committee. Embryos were raised in seawater tanks at 18 C and euthanized by overdose of MS-222 (Sigma) at selected stages (according to Ballard et al. 1993). Zebrafish D. rerio embryos were collected at SciLifeLab Genome Engineering Zebrafish Facility in Uppsala, Sweden, which holds an ethical permit from the Swedish Board of Agriculture. The Ethics Committee of the University of Concepcion approved the experimental procedures used to euthanize and process X. tropicalis tadpoles, which were strictly performed following the recommendations of the Biosafety and Bioethics Manual of the National Commission of Scientific and Technological Research (CONICYT, Chilean Government).

RNA Sequencing
A section of the abdomen including the vertebral column was dissected into RNAlater from one catshark embryo at stage 32 (3.5-month old). The sample was homogenized in liquid nitrogen and total RNA extracted with RNeasy fibrous tissue Mini kit (Qiagen). Total RNA quality and quantity were determined using Bioanalyzer 2100 and Qubit 3.0. The library was prepared from 360 ng total RNA using the TruSeq stranded mRNA library preparation kit, from which the Poly(A)þ mRNA fraction was isolated using poly-T oligo-attached magnetic beads, then chemically fragmented, and finally reverse transcribed to generate cDNAs that are ligated to adapters. These products were purified and amplified with PCR to generate the indexed final double stranded cDNA library. The quality of the library was evaluated in TapeStation 4200, High Sensitivity assay; the quantity of the library was analyzed by real-time PCR in LightCycler 480 (Roche). Three hundred and sixty nanograms of the cDNA library was sequenced by paired-end sequencing (100 Â 2) in Illumina HiSeq 2500 sequencer.
Transcriptome sequencing data were deposited into NCBI SRA database under accession number SRR8753342.

Identification of Collagen Type X Genes in the Genomes and Transcriptomes of Chondrichthyans
Collagen 10a1 (Col10a1), Collagen 8a1 (Col8a1), and Collagen 8a2 (Col8a2) sequences for human, mouse, and zebrafish were collected from GenBank. Elephant shark whole genome assembly GCA_000165045.2 Callorhinchus_milii-6.1.3 (Venkatesh et al. 2014) was downloaded and screened for Col10a1, Col8a1, and Col8a2 sequences using BLAST. Scaffolds with best matching alignments were downloaded and genes were predicted with Genscan (Burge and Karlin 1997) and Augustus (Stanke and Morgenstern 2005). Predicted Col10a1, Col8a1, and Col8a2 protein sequences from elephant shark were used as TBLASTN queries in little skate (L. erinacea) Transcriptomic Contigs-Build 2 of Skatebase (Wyffels et al. 2014) and against the locally assembled lesser spotted catshark and thornback ray jaw and vertebrae transcriptomes on the Montpellier Bioinformatics Biodiversity (MBB) platform (http://mbb.univ-montp2.fr/ MBB//index.php). Nucleotide sequence data reported are available in the Third Party Annotation Section of the DDBJ/ENA/GenBank databases under the accession numbers TPA: BK010855-BK010872. Catshark Col10a1 cDNA sequences were used to search for similar sequences in a whale shark genome (R. typus) using BLASTN (GCA_001642345.2 ASM164234v2 [Read et al. 2017]). Thornback ray and little skate cDNA sequences were used to map synteny on the Thorny skate A. radiata draft genome using BLAST (data accessed with permission from Vertebrate Genome Project https://vgp.github.io/genomeark/Amblyraja_radiata/).

Phylogenetic Analysis
Predicted CDS and protein sequences for all identified collagen type X and type VIII genes from elephant shark, catshark, little skate, thornback ray, and whale shark together with sequences from osteichthyan species were used for phylogenetic tree reconstruction (supplementary table 1, Supplementary Material online). Sequences were aligned with Macse (Ranwez et al. 2011) in order to account for the codon structure, producing separate output files of aligned nucleotide and their corresponding amino-acid sequences. These alignments were then cleaned using hmmcleaner to mask unreliably aligned residues (https://metacpan.org/pod/ HmmCleaner.pl). Phylogenetic analyses were performed on both nucleotide and amino-acid alignments to infer the evolutionary history of these gene families and provide orthology assignments, by rooting the tree with the Col8a1 clade (Aldea et al. 2013). This data set was used to reconstruct gene phylogenies in Maximum Likelihood using IQ-TREE 1.5.4 (Nguyen et al. 2015) under the GTRþR4þF evolution model for nucleotide data and under the LG4XþF evolution model for amino-acid data. Node support was estimated by performing a hundred bootstrap replicates for each analysis.

Protein Domain Prediction
Conserved protein domains were predicted in the sequences of human, elephant shark, and catshark Col10a1 proteins using InterPro (Finn et al. 2017) and illustrated with Adobe Illustrator.

Conserved Synteny Analysis
The synteny analysis was performed on the genomic region including the seven first exons of the Nt5dc1 gene, and on the first exon of the Frk gene, which is the closest coding sequence to the exon 1 of Nt5dc1 in tetrapods ( fig. 1) . The X. tropicalis 9.1 genomic data were visualized from XenBase (Karimi et al. 2018) where the Col10a1 sequence mapped 3 0 to the fifth exon of a gene named Abhd10 which BLAST best hit is Nt5dc1 in the mouse transcriptomic data. Homology between human exon 6 and Lepisosteus exon 9 was verified by CDS alignment. All sequence reference numbers are given in supplementary table 1 (Supplementary Material online). Genomic scaffolds containing identified Col10a1 genes were obtained from genome databases of C. milii (Venkatesh et al. 2014) and R. typus (Read et al. 2017). The six individual Col10a1 cDNA sequences from the catshark transcriptome were used to screen the latest genomic assembly of the catshark using BLASTN on the Galaxy platform on the MBB server (Oulion et al. 2010) which allowed to retrieve six contig sequences including one or more of the Col10a1 cDNA sequences. Catshark Frk and Nt5dc1 coding sequences were identified with TBLASTN searches of the homologous sequences from the elephant shark genome against our transcriptome (references in supplementary table 1, Supplementary Material online). Thornback ray and little skate Col10a1 cDNA sequences were also used to identify synteny on the draft genome of the thorny skate A. radiata using BLASTN.
Sequences were amplified from cDNA extracted from the jaw of stage 33 embryos. PCR products were ligated into the pGEM-T easy vector using the TA cloning kit (Promega). Inserts were amplified using primers present in the vector backbone, and PCR products were used as a template to synthesize antisense DIG riboprobes in a 3 ml reaction containing 100-200 ng PCR product and using the DIG RNA labeling mix (Roche) and T7 or SP6 RNA polymerase (Promega) following manufacturer's instructions. DIGlabelled riboprobes were purified on MicroSpin G50 column (GE Healthcare).
Whole embryos from stage 33 were fixed for 48 h in 4% PFA in 1Â PBS at 4 C and then transferred in 100% ethanol at À20 C for storage. Catshark in situ hybridizations were performed on 14 mm thick cryosections (para-sagittal in lower jaws and transversal in the abdomen) and stored at À20 C. In situ hybridization on sections was performed as described previously (Enault et al. 2015) with stringent conditions of hybridization at 70 C. Images of the catshark in situ hybridization and histological staining were taken with Hamamatsu NanoZoomer 2.0-HT Slide Scanner (Montpellier RIO Imaging facility, INM Optique) with a 40Â objective.
In Situ Hybridization on Zebrafish and Xenopus Larvae Primers were designed for zebrafish col10a1a Fw CCTGGAGCCAAAGGAGAGTT; Rev TATCGGCAGCAA AGACACCA; and zebrafish col10a1b Fw TTCATCTCC TGGGAAGCCTG; Rev TTCACCTCTGCTACCTGGTG gene sequences. Zebrafish cDNA was synthesized with SuperScript III First-Strand Synthesis Kit (ThermoFisher) and oligo dT primers using RNA isolated from zebrafish embryos at 5 days post fertilization (dpf) as a template. Fragments were PCR amplified and cloned from zebrafish cDNA into TOPO pCR II TOPO vector and antisense RNA probes were synthesized with either SP6 or T7 RNA polymerase and digoxigenin labeling mix (Roche). In situ hybridization on zebrafish wholemount larvae was done as previously described (Filipek-G orniok et al. 2013). Whole-mount larvae after in situ were embedded in paraffin and sectioned at 6 mm. Sections were imaged with 40Â objective on Leica DM5500B microscope and Hamamatsu NanoZoomer S60 Digital Slide Scanner.
In situ hybridization on paraffin sections of X. tropicalis stage NF60 tadpoles were performed using previously reported Col10a1 probe and previously described protocol (Aldea et al. 2013