Transcriptomic Analysis Identifies Candidate Genes for Differential Expression during Xenopus laevis Inner Ear Development

Background The genes involved in inner ear development and maintenance of the adult organ have yet to be fully characterized. Previous genetic analysis has emphasized the early development that gives rise to the otic vesicle. This study aimed to bridge the knowledge gap and identify candidate genes that are expressed as the auditory and vestibular sensory organs continue to grow and develop until the systems reach postmetamorphic maturity. Methods Affymetrix microarrays were used to assess inner ear transcriptome profiles from three Xenopus laevis developmental ages where all eight endorgans comprise mechanosensory hair cells: larval stages 50 and 56, and the post-metamorphic juvenile. Pairwise comparisons were made between the three developmental stages and the resulting differentially expressed X. laevis Probe Set IDs (Xl-PSIDs) were assigned to four groups based on differential expression patterns. DAVID analysis was undertaken to impart functional annotation to the differentially regulated Xl-PSIDs. Results Analysis identified 1510 candidate genes for differential gene expression in one or more pairwise comparison. Annotated genes not previously associated with inner ear development emerged from this analysis, as well as annotated genes with established inner ear function, such as oncomodulin, neurod1, and sp7. Notably, 36% of differentially expressed Xl-PSIDs were unannotated. Conclusions. Results draw attention to the complex gene regulatory patterns that characterize Xenopus inner ear development, and underscore the need for improved annotation of the X. laevis genome. Outcomes can be utilized to select candidate inner ear genes for functional analysis, and to promote Xenopus as a model organism for biomedical studies of hearing and balance.


Introduction
The coordinated activity of multiple genes underlies many complex biological processes including organ development, cell differentiation, and the onset and persistence of diseases.High-throughput methods, such as microarrays and RNA-Seq, can quantify the expression levels of thousands of transcripts simultaneously, permitting a system-level analysis of transcriptomes.The use of global profiling technologies enables the discovery of expression profiles that correlate with biological processes, which can result in the generation of hypotheses regarding the genes involved in the maintenance of the phenotype in question (McDermott et al. 2007;Friis et al. 2011;Darville and Sokolowski, 2013;Gu et al., 2016;Baxi et al, 2023;Tisi et al, 2023).
The inner ear is an example of an intricate organ system whose development requires the controlled expression of manifold genes, many of which have yet to be fully defined (Alsina et al, 2009;Wu et al., 2012;Mackowetzky et al., 2021).The inner ear's role as a sensory receptor for mechanical forces imparts its vital function as an organ essential for vertebrate survival and reproduction.The start of the current millennium ushered in the genomics era and subsequently, various genes that are essential for inner ear development and function have emerged.For example, evidence suggests that a number of transcription factors such as DLX5, SIX1, and GATA2 are associated with hearing loss and/or vestibular dysfunction and many Online Mendelian Inheritance in Man ® genes have been linked to deafness and vestibular disorders (Merlo et al., 2002;Zheng et al., 2003;Haugas et al., 2010;Ramirez-Gordillo et al. . 2015;Taiber et al, 2022).However, in comparison with other organ systems such as the eye, heart, and kidney, far less is known about the inner ear's developmental transcriptome, and the control of genetic networks that typify inner ear development (Giraldez and Fritzsch, 2007;Broto et al., 2021;Mackowetzky et al, 2021).Investigations that seek to understand how the expression profiles of functionally relevant genes change during inner ear development and maturation, are necessitated by the pressing need to develop therapeutic interventions for disorders of hearing and balance that affect almost half a billion humans (Agrawal et al., 2009;Eshraghi et al., 2013;Lustig and Akil, 2019;Tabier et al., 2022).It is anticipated that transcriptomic profiling will uncover crucial inner ear genes whose function is relevant for prevention or repair of malfunctioning auditory and vestibular systems, as well as for in vitro studies of inner ear organs in stem cell derived systems ( Koehler and Hashino, 2014;Lee and Waldhaus, 2022).
The inner ear's unique role, sequestered location within the temporal bone, and relatively minute size make studies with human tissue impractical.Consequently, investigators who study hearing and balance rely on a spectrum of animal models and advocate for the advantages of their particular organism as a model for investigating inner ear structure and function (Giraldez and Fritzsch, 2007;Mackowetzky et al, 2021;Lee and Waldhaus, 2022).
Common justifications for species selection include the organism's relevance as a mammalian model (rodents, cats, primates), its regenerative potential (fish, amphibians, birds), or its suitability for transgenesis (Danio Rerio, Xenopus, Mus musculus).In fact, foundational understanding of the physiological role of the inner ear's mechanoreceptor hair cells was pioneered in non-mammalian models through seminal biophysical investigations of amphibian and turtle hair cells (Hudspeth et al., 1977;Ricci et al., 2002;Kozlov et al. 2007).
Research presented here implements the amphibian, Xenopus, for developmental investigations of inner ear transcriptomics using microarray approaches.A popular and important genus for studies of embryogenesis (Gurdon, 2014;Kostiuk and Khoka, 2021), Xenopus has been relatively underutilized for studies of organogenesis, neural systems, or aging, especially in comparison with zebrafish, an organism that is evolutionarily far more distant from humans than amphibians (Nakatani et al., 2007;).With two sequenced genomes (X.laevis and X. tropicalis),

3
Resource Centers that can produce and stock transgenic animals (Vize and Zorn, 2017;Horb et al, 2019) and documented neuroregenerative potential (Burns et al., 2013;Phipps et al., 2020), Xenopus is poised to contribute to genomic investigations of neural and sensory processes (Nenni et al., 2019;Exner and Willsey, 2021) .Previous transcriptomics studies featuring Xenopus have profiled multiple organs and developmental stages, but few have focused on the inner ear and its development (Baldessari et al., 2005;Powers et al., 2012;Langlois and Martyniuk 2013).The paucity of global profiling studies that feature the inner ear represents a missed opportunity to identify new genes for inner ear and organ function.Technologies that measure transcript abundance across multiple genes also afford the opportunity to discover groups of genes with synchronized changes in expression.Identification of such genes has the potential to uncover previously unknown functional relationships.
We implemented a microarray approach for analysis of the Xenopus inner ear transcriptome at three developmental ages that represent anatomically distinct stages of inner ear organ formation: s50-s52 (s50), s56-s58 (s56), and juvenile (aged 3 months).Larval age s50 represents the youngest developmental age when all auditory and vestibular end organs have been anatomically characterized as fully formed (s50) (Bever et al., 2003;Quick et al., 2005).The s56 age was selected as an intermediate larval age during which the inner ear undergoes an expansive growth period as evidenced by both the development of new hair cells and innervation by the eight cranial nerve (Díaz et al., 1995;Lopez-Anaya et al., 1997;Bever et al., 2003).Finally, the juvenile age was selected to identify alterations in gene expression related to the structural and functional changes that occur after metamorphosis, when Xenopus enters reproductive age.For example, although animal size increases after metamorphosis, very few hair cells are added to the sensory field and the 8 th nerve undergoes extensive myelination as compared with larval ages (Diaz et al, 1995;Lopez-Anaya et al., 1997).
Differential expression analysis of pairwise comparisons between developmental stages enabled identification of transcripts with the most variable expression levels between ages, as well as those that were predominantly expressed during a single age.The genes identified from the pairwise comparisons were subjected to functional annotation clustering using the Database for Annotation, Visualization and Integrated Discovery (DAVID); transcripts upregulated in the juvenile in comparison to the two larval stages were categorized as the cohort with the greatest enrichment score.The analysis presented here emphasized identifying the genes with the greatest expression changes to produce a robust data set that expands the characterization of three distinct developmental stages of the Xenopus inner ear.Research outcomes enhance the annotation of the Xenopus genome, and offer a resource that can be used by the inner ear development community to generate hypotheses about sensory organ formation, as well as for investigators interested in drug discovery (Wheeler and Liu, 2012;Liu and Yang, 2022;Cousins, 2022), and the use of stem cells for investigations of inner ear organ systems (Zhang and Hu, 2012;Roccio, 2021).

Evaluation of replicate data for technical consistency of sample and array methods
Raw, log2 transformed, data for each of the three developmental stages and the associated replicates are shown in the boxplots in Fig. 1A.Visual inspection of boxplots prior to preprocessing showed similar medians in stage replicates, with the s56 replicates appearing the least variable.In order to achieve the broadest intensity distribution, raw array data were preprocessed using GCRMA.Preprocessing resulted in normalized data with reduced variability in comparison to raw data (Fig. 1B).Examination of scatter plots of all samples plotted against 4 one another (Fig. 1C), showed the greatest correlation to be between replicate samples and the least between different developmental stages, with the exception of the s50 replicates.Histograms of all developmental stages were also constructed.The histograms displayed similar intensity distributions with the greatest intensity peak at the lowest observed intensity level (Fig. 1C), demonstrating the sensitivity of the GCRMA method to low intensity levels in comparison to other preprocessing methods such as RMA (Wu et al., 2004).It was concluded that the stringent surgical procedures and RNA isolation methods utilized in this study (Trujillo-Provencio et al., 2009) enhanced reproducibility and minimized technical experimental variation, based on (1) the similarities observed in both the replicate raw and preprocessed data, and (2) the lower observed variability between replicates than between developmental stages.

Differential expression analysis identified 1510 differentially expressed Xl-PSIDs (DEX) that were assigned to seven differential expression categories
The objective of this study was to determine gene expression changes during inner ear development using Xenopus laevis as a model.To achieve this objective, pairwise comparisons between the three developmental stages were made using the convention that a positive fold change specifies greater expression in the older stage.Data that met the filter criteria (q-value ≤ 0.01; fold change ≥ 1.5) were used for analysis.The filter criteria were implemented to focus on differences with a greater likelihood of biological significance, but to also include changes between genes that may be represented at low abundance such as transcription factors.
Volcano plots of the comparisons are shown in Fig. 2 A-C, with the dotted line crossing the x-axis at the significance level of q-value ≤ 0.01.The stages analyzed in this study included two different larval stages (s50 and s56) and the metamorphosed juvenile (Juv).The fewest differentially expressed Xl-PSIDs (DEX) were found in the pairwise comparison featuring the two larval stages (Fig. 2A).The proportion of significant Xl-PSIDs increased in the Juv-s56 comparison (Fig. 2B) and, finally in the Juv-s50 plots (Fig. 2C), more data points were significant than in any other comparison.The heatmap shown in Fig. 3A displays hierarchical clustering and expression differences between all samples after filtering by statistical significance and fold change.In the heatmap, all replicates clustered by stage and the two larval stages clustered together.
The resulting DEX could be assigned to seven categories based on their expression pattern (Fig. 3B).The division into the seven differential expression categories served as the overall focal point of the downstream analysis.The most populated categories comprised Xl-PSIDs that were differentially expressed in one (Group A) or two (Group B) pairwise comparisons, these two groups accounted for 98% of the DEX.The least populated category (Group C) encompassed Xl-PSIDs differentially expressed between all three stages, with 2% of the DEX total.

Fig. 2. Results of pairwise comparisons between the ages. (A,B,C) Volcano plots of stage pairwise comparisons:
s56-s50 (A), Juv-s56 (B), and Juv-s50 (C).The Juv-s50 plot displays the greatest variation in the magnitudes of the log2 fold change in the differential analysis.Dotted line represents -log10(p) = 3.12; data points above this line are statistically significant at q ≤ 0.01.with the greatest contribution to the total arising from Xl-PSIDs differentially expressed in the both the Juv-s50 and Juv-s56 comparisons.The threshold for differential expression was set at q-value ≤ 0.01 and minimum fold change difference of 1.5 a.u.

Characteristics of Xl-PSIDs differentially expressed in a single pairwise comparison (Group A)
The greatest number of Xl-PSIDs in Group A were noted in the Juv-s50 category, followed by the Juv-s56 and s56-s50 comparisons (Fig. 3B).In Group A, 34-38% of Xl-PSIDs differentially expressed in a single comparison were unannotated.The majority (61%) of the 159 DEX in the s56-s50 comparison were upregulated in s50 (UP_s50) and the observed fold changes ranged from -4.6 to 6.4 (Fig. 4A, Table 1).Among the most upregulated Xl-PSIDs in this comparison, two were unannotated and three were annotated with the following gene symbols: atp6v1g3, LOC398308 (lectin type 2), and tmod4.Two Xl-PSIDs among the top expressed in s56 also lacked annotation.The remaining four from the top UP_s56 Xl-PSIDs corresponded to the gene symbols: cldn1, slc25a13, klf9-b, and hhip.The Xl-PSID annotated as cldn1 (claudin 1), with a fold change of 4.9, was the Xl-PSID with the second highest fold change in this group (the highest fold change was for an Xl-PSID that was not annotated).The hedgehog interacting protein (hhip) is an inhibitor of hedgehog signaling, which is involved in anterioposterior patterning in Xenopus ( Waldman et al., 2007).
In the Juv-s56 comparison, in contrast to the s56-s50 comparison, greater upregulation (67%) was seen in the more mature developmental stage (Fig 4B).This comparison contained the second greatest count of DEX seen in this study and the fold change range was greater than in s56-s50, ranging from -10.7 to 9.6 (Table 2).The two Xl-PSIDs UP_56 with the greatest fold changes (-10.7 and -8) were annotated as mediator complex subunit 22 (med22) and major histocompatibility complex, class I, A (hla-a).The Xl-PSID annotated as deiodinase, iodothyronine, type 2 (dio2) had the fourth highest fold change of -5.8 and is involved in cochlea development in mice by way of local regulation of thyroid hormone (Campos-Barros et al., 2000).The Xl-PSID UP_Juv with the greatest fold change was identified as calpastatin-like (calp1), which is an inhibitor of a family of calcium activated proteases known as the calpains, which are involved in apoptotic processes (Rojas et al., 1999).
In contrast to both the s56-s50 and Juv-s56 comparisons, the Juv-s50 comparison had a similar proportion of upregulated Xl-PSIDs between the older (54%) and younger ages (Fig 4C).In this comparison, the Xl-PSID foldchanges ranged from -9.8 to 11.9 (Table 3).The Xl-PSID most upregulated in s50 was annotated as pavlb.2and the second most upregulated was oncomodulin (ocm.2), with a fold change of -9.1.Oncomodulin, also known as parvalbumin-β, is a member of the parvalbumin family of calcium binding proteins that are expressed in the mammalian cochlea and differentiates outer from inner hair cells (Sakaguchi et al., 1998).The expression level of oncomodulin decreases as development progresses in the rat inner ear (Yang et al., 2004), which was also observed in this analysis.One of the UP_s50 Xl-PSIDs was transcription factor neurod1-b, with fold change of -6.6.The transcription factor neurogenic differentiation 1 (neurod1) is an important factor in the development of both the cochlea and the vestibular apparatus (Liu et al., 2000) and it represses the transformation of neuronal sensory cells into hair cells (Jahan et al., 2010).The most upregulated Xl-PSIDs in the juvenile stage in this comparison were the hemoglobins hbg2-a and hbg1, with fold changes of 11.9 and 10.3 respectively.

Analysis of Xl-PSIDs differentially expressed in two comparisons facilitated a "stage_centric" view of differential expression (Group B)
Xl-PSIDs that were differentially expressed in two comparisons (Group B) were termed stage_centric.
Using the s50 age as an example, the set of Xl-PSIDs differentially expressed when s50 was compared to both s56 and juvenile would represent the group of s50_centric Xl-PSIDs.
Analysis of DEX in s56-s50 ∩ Juv-s50 (s50_centric) provided information regarding differentially regulated genes in s50 relative to a more developed larval stage and the completely metamorphosed juvenile.Genes that are UP_s50 in this context could represent a cohort with decreasing expression during development.Most (55%) of the Xl-PSIDs in this category were upregulated in s50 (Table 4).There were 215 DEX in Juv-s50 ∩ s56-s50, 36% of which did not have a gene annotation.The Xl-PSID with the greatest upregulation in s50 was annotated as dentin sialophosphoprotein (dspp), with fold changes of -8.9 (Juv-s50) and -6.3 (s56-s50) [Table 5].Mutations in dspp have been implicated in dentinogenesis imperfect, which can present with an autosomal dominant form of hearing loss (Xiao et al., 2001).Also among the most upregulated s50 Xl-PSIDs were several genes with functional roles in cellular structure, two of which were the keratin genes krt.5.6 and krt14, both with log2 fold changes of -8.
The annotated Xl-PSID with the greatest downregulation in s50 was annotated as bone gamma-carboxyglutamate [gla] protein (bglap), with fold changes of 8.5 and 7.8, compared to juvenile and s56, respectively.
Xl-PSIDs that were differentially expressed in s56-s50 ∩ Juv-s56 comparisons (s56_centric) comprised the lowest count of DEX in Group B comparisons with a total of 98 (Fig. 3B), 49% of which were upregulated in s56.
The s56_centric cohort of genes may represent differentially expressed genes whose expression peaks as metamorphosis approaches, and then decrease after its completion.This gene set could include genes with a role in preparing the inner ear for metamorphosis.In addition, the percentage of Xl-PSIDs without annotation was greater in this category than in s56_centric at 40%.Two of the top 15 Xl-PSIDs upregulated in s56 were transcription factors (Table 6).The transcription factor annotated as specificity protein transcription factor (sp7) was among the most upregulated in s56 with fold changes of -5.6 (Juv-s56) and 7.8 (s56-s50).A genetic mutation in sp7 has been linked to osteogenesis imperfecta (Lapunzina et al., 2010).The other transcription factor in this category was the Xl-PSID annotated as distal-less homeobox 3 (dlx3-b), which was one of the lower expressed of the top 15 with fold changes of -4 (Juv-s56) and 3.2 (s56-s50).The three Xl-PSIDs with the greatest downregulation in s56 were unannotated.The most downregulated annotated Xl-PSID was calcium binding and coiled-coil domain 1 (calcoco1) with fold changes of 4.9 (Juv-s56) and -3.4 (s56-s50).The calcoco1 gene is considered a positive regulator of transcription and has been shown to be expressed in multiple tissues including the brain and kidney (Kim et al., 2003).
The juvenile centric Xl-PSIDs, as shown in the Venn diagram in Fig. 3B, represented the most populated category with 432 Xl-PSIDs.This category had the fewest unannotated Xl-PSIDs at 33%, and 48% of the Xl-PSIDs were upregulated in the juvenile (Table 4).The two most upregulated Xl-PSIDs in juvenile were annotated as the same hemoglobin gene (hba1), with fold changes ranging from 8.7 to 11 (Table 7).Another Xl-PSID highly upregulated in juvenile was annotated as solute carrier family 25 [mitochondrial carrier; adenine nucleotide translocator], member 5 (slc25a5), with a fold change of 9 in both Juv-s50 and Juv-s56.
Hemoglobin genes were also among the Xl-PSIDs with the greatest downregulation in juvenile.Twelve of the top 15 Xl-PSIDs downregulated in juvenile were annotated, six of which were annotated as hemoglobin genes.In addition, three of the annotated Xl-PSIDs were annotated as the same hemoglobin gene (hbg2-b), with fold changes ranging from -9.2 to -11.1.

Analysis of Xl-PSIDs differentially expressed in all three comparisons identified genes with changing expression profiles throughout development (Group C)
This expression category contained the fewest DEX and contained the greatest proportion of unannotated Xl-PSIDs at 55%.There were equal percentages (41%) of DEX in all three stages with expression profiles that decreased or increased during development (Fig. 4 D-F).The annotated Xl-PSIDs that decreased during development were hemoglobin, epsilon 1 (hbe1); keratin 5, gene 6 (krt5.6);keratin 12 (krt12); aurora kinase b (aurkb-b); and protein kinase, cGMP-dependent, type II (prkg2) [Table 8].The cohort of annotated Xl-PSIDs with increasing expression during development included dipeptidyl-peptidase 4 (dpp4); Corticotropin releasing hormone receptor 1, gene 2 (crhr1.2);and CD302 molecule (cd302).The other expression pattern observed had intensity values that increased from s50 to s56, then decreased from s56 to juvenile; 18% of the Xl-PSIDs possessed this pattern and the only one with a gene annotation was the presumed transcription factor Kruppel-like factor 5 (klf5).The gene klf5 was also among the most differentially regulated genes in Group B, but was associated with a different Xl-PSID.

Xl-PSIDs for genes encoding proteins with confirmed inner ear function were identified as differentially expressed
High-throughput technologies are powerful tools for transcriptomic analysis because they enable many genes to be profiled simultaneously, thereby facilitating the analysis of the relationships between genes with similar expression profiles.High-throughput technologies also provide a global perspective of gene expression in the system under investigation that can offer a rich context for focused studies of individual or smaller cadres of genes (Baldessari et al., 2005;McDermott et al. 2007;Darville and Sokolowski, 2013).The quality of the data is an essential component of experimental design that must be considered when generating and analyzing high-throughput data.Accordingly, we implemented a protocol to produce a microarray data set that was optimized to enhance reproducibility and minimize experimental variation generated by technical procedures.Evidence of this optimization can be seen in Fig. 1A, where the similarities in intensity distribution of the replicate data are apparent even prior to normalization.
Our microarray data analysis identified 1510 distinct DEX.The DEX included genes with established function in the inner ear, comprising genes associated with hearing loss and/or vestibular dysfunction such as transcription factors sp7, and neurod1-b.Neurod1 is expressed during otic vesicle specification and is known to be involved in both auditory and vestibular development (Liu et al., 2000;Alsina et al., 2009).Mutations in neurod1 are linked to a clinical syndrome whose symptoms include sensorineural deafness (Rubio-Cabezas et al., 2010).In addition, genes were identified that are related to syndromic hearing loss such as structural proteins col1a1 and col1a2 (implicated in Osteogenesis imperfect [Marini et al., 2007]) and the transcription factor tfap2a-b (linked to Branchio-oculo-facial syndrome [Tekin et al., 2009]).The detection of known inner ear genes in this data set enhanced confidence in the technique and demonstrated that transcripts detected in the inner ears of other vertebrates were also detected in the Xenopus inner ear.

Pairwise comparisons uncovered Xl-PSIDs with seven differential expression patterns, including those for genes not previously associated with inner ear function
Pairwise comparisons between the three developmental stages produced seven differential expression patterns that were partitioned into Group A, B, or C contingent on whether an Xl-PSID was differentially expressed in one, two or three pairwise comparisons; respectively.Of these Xl-PSIDs, 64% corresponded to genes with an established functional annotation.Analysis of Group A Xl-PSIDs showed upregulation of gene expression in the s50 and in the juvenile stages, and downregulation of expression in s56, with the greatest differences observed between s50 and juvenile stages (Fig. 4 A-C).Although the inner ears of the s50 and juvenile animals are both functional, dramatic changes occur during this period such as an increase in the number of hair cells, axon projections, and in the size of the inner ear overall that are all likely contributors to the differences in gene expression observed (Díaz et al., 1995;López-Anaya et al., 1997;Serrano et al., 2001).
A concept implemented in this study is that of "stage centricity", which designates sets of Xl-PSIDs in which the direction of expression in one stage (the centric stage) is in opposition to that of the other two stages.The comparisons that produced stage centric expression patterns were designated as Group B (Table 7, Fig. 3B).In Group B, the s50_centric Xl-PSIDs had the greatest proportion of upregulation with 55% upregulated in s50, indicating greater gene expression in the youngest stage when compared to both the s56 and the juvenile.The juvenile_centric expression category comprised 29% of the Group B Xl-PSIDs, twice as many as s50_centric and four times the number in s56_centric.This further supports the conclusion that a larger difference in gene expression is found between the juvenile and larval stages than between the two larval stages.In addition, there appears to be more downregulation of gene expression as the inner ear ages from larval to the juvenile stage.
The smallest cohort of Xl-PSIDs were found in Group C (Fig. 4 D-F, Table 8), which comprised the 22 Xl-PSIDs that were differentially expressed between all pairwise comparisons.Twelve of the Group C Xl-PSIDs did not have a gene annotation, making the most actively regulated Xl-PSIDs the least annotated category in this analysis.The dynamic nature of the genes that are changing expression throughout the developmental stages examined combined with the lack of annotation make these Xl-PSIDs of particular interest since their fluctuating expression levels may indicate a functional role in inner ear development and maturation that has yet to be characterized.

Xenopus gene annotation limits interpretation of differential expression patterns
DAVID functional annotation clustering was undertaken on Xl-PSIDs from Groups A & B as summarized in Table 9.A total of seven significant clusters were discovered, five of which were found for Xl-PSIDs downregulated in juvenile in the juvenile_centric expression category.We noted that Xl-PSIDs downregulated in juvenile in the juvenile_centric set were both the most abundant (224 Xl-PSIDs), and had the greatest proportion of Xl-PSIDs that mapped to DAVID IDs (57%).A recurring theme encountered in this study was a lack of annotation for many of the differentially regulated genes.Three of the seven differential expression categories that produced significant clusters included terms previously associated with inner ear function or development, but the inability to map many of the Xl-PSIDs to DAVID IDs likely had a negative impact on the effectiveness of the functional analysis.

Toward a Xenopus model for hearing and balance
The inner ear remains a relatively understudied organ.Data from this Xenopus transcriptome-wide study can be mined to identify new candidate genes for inner ear function and to explore unrecognized relationships between inner ear genes.The dataset presented here complements and extends genetic outcomes from inner ear research that has focused on induction and formation of the otic placode (Giraldez and Fritzsch, 2007;Alsina et al, 2009;Almasoudi and Schlosser, 2021)..The observation that approximately a third of the differentially expressed Xenopus inner ear genes have no annotation highlights the prevalence of knowledge gaps regarding genes that are involved in sensory organ development.These unannotated genes represent an exciting opportunity to fully characterize the genetics underlying inner ear development and to identify genes crucial for balance and audition.
Additionally, when genetic characterizations are carried out in an amphibian like Xenopus, we have the future prospect of identifying genes involved in hair cell regeneration, a process that does not occur easily in the mammalian inner ear (Costa et al., 2015;Lee and Waldhaus, 2022).It is therefore anticipated that the data set presented here will become a resource that can be used to generate hypotheses regarding genes and mechanisms that underlie inner ear development and repair.
In summary, Xenopus laevis is an established model for vertebrate development and cellular biology, especially during embryogenesis (reviewed by Harland and Grainger, 2011;Gurdon, 2014), but its allotetraploid genome historically was not easily amenable to genetic manipulation.Previously this limitation was overcome using

Preprocessing for microarray Analysis
The original (raw) data from X. laevis GeneChip ® CEL files acquired from three replicate arrays per developmental stage were preprocessed using the Gene Chip robust multichip averaging (GCRMA) method to produce a single log2 transformed measure of the intensity level (in arbitrary units [a.u.] of fluorescence) for every Xl-PSID on each replicate array.The JMP Genomics 5.0 analysis platform was used for both preprocessing and differential expression analysis.The preprocessed microarray data served as the starting point for all downstream analysis.Data can be accessed at the Gene Expression Omnibus under accession numbers GSE69546, GSE73828, GSE73829.

Differential expression analysis and pairwise comparisons
Pairwise comparisons were made between all stages using Analysis of Variance (ANOVA) as implemented by the JMP Genomics 5.0 analysis platform.Candidate Xl-PSIDs for differential expression were detected by analyzing pairwise differences in the average fluorescence intensity of the replicate Xl-PSIDs for each stage.
Filters and fold change restrictions were applied to the preprocessed data to identify candidate genes for differential expression during development.A q-value (≤ 0.01) established criteria for a positive false discovery rate (pFDR).The minimum fold change value for pairwise comparison was set to ±1.5 to allow for the capture of subtle changes in expression and to filter out potentially insignificant changes.As per Powers et al. (2012), only Xl-PSIDs that met the additional minimum average normalized intensity criteria of fluorescence intensity ≥ 4 A.U. on the replicate arrays were designated "upregulated".
The resulting data were exported as a tab-delimited text file for the application of significance filters.A negative fold change in a given pairwise comparison was used to identify downregulated Xl-PSIDs.Java software was developed to apply the filters against the JMP output file.Following application of the fold change and fluorescence intensity significance filters, the differentially expressed data were stored in a MySQL database table.
The Java software developed for this project will be made available by the authors and deposited to the GitHub open source project repository (https://github.com/).

Annotation of Xl-PSIDs
Xl-PSIDs were annotated using a two-stage approach.Data files were downloaded from the Xenbase

Fig. 1 .
Fig. 1.Comparison of microarray data distributions for samples of different ages.(A, B, C) Box plots of Affymetrix Genome 2.0 array data before (A) and after preprocessing with GCRMA (B).(C) Scatterplots of the Xl-PSID intensity values from the three pairwise comparisons show greatest similarity between replicate comparisons of the same animal age.Fluorescence intensity values are in arbitrary units of fluorescence (a.u.).

Fig 3 .
Fig 3. Heatmap and Venn Diagram of the 1510 differentially expressed Xl-PSIDs.(A) The heatmap shows the results of hierarchical clustering of all developmental stage replicates.Clear differences are visible between the developmental stages and the replicates from the two larval stages are clustered together.(B) Venn diagram depicting the distribution of Xl-PSIDs differentially expressed in one, two, or all three stage pairwise comparisons,

Fig. 4 .
Fig. 4. The smallest differential expression cohort was comprised of the 22 Xl-PSIDs differentially expressed in all three pairwise comparisons.(A-C) Microarray analysis uncovered 743 Xl-PSIDs that met the criteria for differential expression in only one pairwise comparison.In comparisons featuring the juvenile stage, the majority of Xl-PSIDs were upregulated in this stage as opposed to either of the larval stages.(D-F) The Xl-PSIDs differentially expressed in all comparisons had greater upregulation in s56 when compared to both the larval s50 and juvenile stages, in contrast to the pattern observed in Xl-PSIDs differentially expressed in only a single comparison.Circle size adjusted to reflect number of genes.
website (www.xenbase.org;Bowes et al., 2010) and formatted for insertion into a local MySQL database to facilitate automated queries.Differentially expressed Xl-PSIDs were linked to Xenbase annotation by searching data originating from the Xenbase data file GenePageAffymetrix_laevis2.0.txt (file date stamped January 29, 2015) which associated Affymetrix Xl-PSIDs to Xenopus gene symbols.Xl-PSIDs that remained unannotated by Xenbase were subjected to a second round of annotation process.The representative public id from the Affymetrix annotation file, X_laevis_2.na33.annot.csv,was queried against the Xenopus laevis UniGene build 94 Xl.data file for retrieval