Teleost exhibit remarkable diversity of sex determination mechanisms, including genetic sex determination (GSD; Myosho et al., 2012), environmental sex determination (ESD; Li et al., 2023), and various combinations of these two modes. In the GSD system (such as XX/XY and ZW/ZZ), the expression of sex-determining genes initiates a series of cascadic signaling pathways of sex determination and differentiation, thereby inducing the primordial gonads to develop into ovaries or testes (Li and Gui, 2018). Many master sex-determining genes or sex differentiation-related genes have been reported to play a main role in regulation of sex development in various fishes, such as dmrt1 (doublesex and mab-3 related transcription factor 1), sox9 (SRY-box transcription factor 9), and foxl2 (forkhead box protein L2). Some transcription factors, involved in regulating gene expression, play a crucial role in the development of male or female differentiation (Cocquet et al., 2002; Cui et al., 2017; Hu et al., 2021).

In fish, dmrt1 is the master sex-determining gene in half-smooth tongue sole (Cui et al., 2017), and its homologous gene dmy (the dmrt1 paralogous located in Y-chromosome) was verified to be the sex-determining gene in medaka (Nanda et al., 2002; Kobayashi et al., 2004). Nearly all polypeptide chains encoded by dmrt genes contain a highly conserved zinc finger DNA-binding and transcription regulating motif (known as DM domain), which consists of six conserved cysteines and two histidines for binding into the minor groove of any target DNA. High expression of dmrt1 in males activates genetic programs to promote testicular differentiation. There is only one copy of the dmrt2 gene in mammals, but two copies (dmrt2a and dmrt2b) are commonly identified in teleost fishes. It has been reported that the dmrt2 (a/b) genes are involved in determination of the body axis and segment differentiation of vertebrates and are also related to the gonadal sex differentiation process (Han et al., 2021). In bony fishes, dmrt3 is usually highly expressed in the testis and nervous system, which presumably plays an important role in the differentiation and development of nerve and germ cells as well as testes (Li et al., 2008; Dong et al., 2010). In contrast, in females, high expression of foxl2 inhibits the expression of dmrt1 while upregulates cyp19a1, thereby promoting the development of ovarian tissues (Zhang et al., 2017). The transforming growth factor-β (TGF-β) signaling pathway mainly participates in mediating the formation of tissues and organs as well as reproductive development by regulating cell growth, proliferation, differentiation and other processes (Li and Gui, 2018; Pan et al., 2021). Many previous studies have proved that this signaling pathway is related to sex determination and differentiation in various fishes. Meanwhile, some members of the TGF-β superfamily have been considered as initiators or key regulators of sexual differentiation, such as amh (anti-Müllerian hormone), amhr2 (anti-Müllerian hormone receptor type 2; (Kamiya et al., 2012), and gsdf (gonadal soma-derived factor; (Myosho et al., 2012).

Transposable elements (TEs) have been integral parts of vertebrate genomes throughout the evolution process and play important roles in promoting genetic innovation, sex determination and reproduction (Biémont and Vieira, 2006). Recent bursts of TEs and simple repeat accumulations were observed around young sex determination loci, indicating the strong association between TE amplification and the formation of sex determination regions on sex chromosomes (Chalopin et al., 2015). TEs could carry regulatory elements and modify the expression of neighboring genes, thereby participating in the evolution of germ cells and gonadal regulatory networks. The expression and transposition of TEs are critical for vertical transmission to progeny and persistence in lineages (Dechaud et al., 2019). Interestingly, TEs are particularly prone to being recruited in sexual development, since they are typically expressed in the gonads (Brunet et al., 2018). For instance, in Oryzias latipes approximately 1.2% and 3.5% of ovarian and testicular transcripts were reported to be attributed to TE expression, respectively (Dechaud et al., 2019).

Asian arowana (Scleropages formosus, also known as dragonfish) belongs to the ancient family Osteoglossidae, and it is a freshwater species with a native distribution in Southeast Asia. It has become a well-known ornamental fish for its vibrant colors and cultural significance. Mainly due to habitat destruction and overfishing, arowanas have been classified as endangered species in Appendix I by the Convention on International Trade in Endangered Species of Wild Fauna and Flora (CITES). In recent years, significant progress has been made in the genetics research and breeding techniques of arowanas (Chang, 2009; Mu et al., 2012; Shen et al., 2014; Austin et al., 2015; Bian et al., 2016; Mu et al., 2022).

Asian arowana often reaches sexual maturity at 2–4 years of age (Chang, 2009), and several previous studies proposed that the genetic sex determination system of this species may be ZW/ZZ (Shen et al., 2014; Bian et al., 2016). Because there are no obvious morphological differences during all life stages of the Asian arowana, even after sexual maturity, it is very difficult to phenotypically differentiate the sex of individuals. One recent report showed that its putative ZW system is more likely a polymorphic pattern that occurs in the 18th chromosome pair, including the accumulation of constitutive heterochromatin and 18S rDNA (Toma et al., 2023). Meanwhile, although its genomes of both sexes have been published by us (Mu et al., 2022), our in-depth comparative genomics between males and females has failed to identified any significant region with obvious genomic differences. Therefore, we performed this study to identify coding and non-coding regions that are differentially expressed between ovaries and testes of adult Asian arowanas through transcriptome sequencing, which may become differential signatures between both sexes for practical gender identification and one-to-one captive breeding.

Transcriptomics is an effective method to obtain gene regulatory networks among different individuals. In this study, we employed the RNA-Seq technology for transcriptomics analysis, and many genes related to gonad development were identified and characterized to improve our understanding of their sex-related functions from a molecular perspective. Generally speaking, the fish sex determination mechanism usually involves complex biological processes, including a series of genes that promote or maintain the development of gonads into testes or ovaries. The sexual dimorphism in expression identified in this study, especially for those previously reported candidate genes that are related to sex determination and sexual differentiation (Kitano et al., 2023), may be effective indicators for sex prediction in Asian arowana.

The Wnt signaling is a crucial pathway to regulate cell proliferation, differentiation, embryonic development, and folliculogenesis. Recent genetic research has identified the wnt4 as a vital regulator of the Wnt/β-catenin signaling pathway during sex determination (Kossack et al., 2018; Farhadi et al., 2021). The most important member wnt4 gene has been identified in zebrafish (Kossack et al., 2018), medaka (LI et al., 2012) and various other teleost fishes. Due to the teleost-specific whole-genome duplication (WGD) events, there are two paralogs in most teleost, wnt4a and wnt4b, both of which play essential roles in the development of female gonads, ovarian growth, early male development, and sex transition (Hu et al., 2014). We found that several wnt genes showed sex-biased expression patterns in arowana. For example, the expression of wnt4a (log2FC = −7.26) and wnt11 (log2FC = −5.46) that are related to ovary differentiation were upregulated in ovaries; the expression of wnt6 (log2FC = 7.95), which could initiate Wnt/β-catenin signaling to control the proliferation of undifferentiated spermatogonia in mouse (Takase and Nusse, 2016), was also testis-biased. rspo1 showed a sexually dimorphic expression pattern with significantly higher expression in testes over ovaries in our present study (log2FC = 5.41). It is a potential female-determining gene in mammals that can regulate the important Wnt/β-catenin signaling pathway; loss of rspo1 can lead to female sex reversal. However, the expression pattern of rspo1 in fish seems not very conservative. For instance, in many teleost fishes such as zebrafish, medaka, and half-smooth tongue, rspo1 shows a female bias, while in some ancient fishes (like Acipenser stephensi, coelacanth, and lungfish) and some sex-changed teleost fishes (like Notolabrus celidotus) shows a male bias(Biscotti et al., 2018; Zhang et al., 2020; Muncaster et al., 2023), implying that rspo1 may have other gender-related functions in various fishes.

The TGF-β signaling pathway is a large family comprising many members responsible for regulation of various processes including tissue and organ formation, reproductive development, cell growth, proliferation, and differentiation (Pan et al., 2021). In fish, this critical pathway often plays a key role in regulating important physiological processes, such as embryonic development, tissue regeneration, immune regulation, chondrogenesis, as well as sex determination and sexual differentiation (Pan et al., 2021; Kitano et al., 2023). Some members belonged to the TGF-β superfamily, such as amh, amhr2, and gsdf were highly expressed in males (Figure 3A), which indicates that these pathways should play a significant role in sexual differentiation of Asian arowana. During ovarian development, the interaction among granulosa cells, theca cells, and oocytes is crucial for follicular development, steroidogenesis, and female oogenesis (Richards, 2018). The forkhead-box L2 (foxl2) transcription factor is one of the earliest known markers of granulosa cell differentiation (Cocquet et al., 2002), and it has been determined to be a major gene of sex determination and maintenance in various teleost (Cui et al., 2017; Zhang et al., 2017; Yuan et al., 2021). In addition, the foxl2 can directly bind to the promoter region of the aromatase gene cyp19a1a to promote its expression, or indirectly regulated expression of cyp19a1a through interaction with sf1 (splicing factor 1), resulting in estrogen (E2) production to promote ovarian differentiation (Wang et al., 2007; Zhang et al., 2017; Yuan et al., 2021). In our present study, folx2 (log2FC = −3.67) and cyp19a1 (log2FC = −6.97) were found to be highly upregulated in the ovary compared with the testis (Figure 3A), showing a similar expression pattern to other fish species (Zhang et al., 2019; Mustapha et al., 2022). Through a PPI analysis, we found that most of the potential sex-related genes among the DEGs exhibit interaction relationships, suggesting that these genes may play functional roles in sex differentiation and development, potentially being directly or indirectly involved in gonadal development or differentiation. Notably, dmrt1, amh, as well as cyp19a1 (and its upstream regulator foxl2) interact with multiple detected potential sex-related genes, indicating that these gene play central regulatory roles in sex differentiation in Asia arowana. In mammals, Dmrt1 directly represses Foxl2 expression in the testis, and Foxl2 is required for repressing Dmrt1 expression in the ovary (Matson et al., 2011). Previous studies on tilapia and zebrafish found that dmrt1 has an antagonistic effect on the expression of foxl2 and cyp19a1 (Wang et al., 2010). Dmrt1 can directly repress the expression of foxl2 and cyp19a1 in somatic cells of tilapia, while in females foxl2 activates the expression of cyp19a1 and inhibits expression of dmrt1. In addition, studies have found that foxl2 has ovary-enriched expression in zebrafish, just like in mammals, and dmrt1 is required for its downregulation. Therefore, it is speculated that the antagonistic relationship between foxl2 and dmrt1 in the sex determination process may be conserved among vertebrates (Webster et al., 2017). We have observed that dmrt1 not only interacts with male-biased gene, such as amh and sox9, but also engages with cyp19a1 and foxl2, suggesting that the antagonism may also exist in Asian arowana. However, due to limited research on Asian arowana, some protein-protein interactions may require further experimental validation (like Dmrt2a). According to the sex-biased DEGs detected in our present study and the sex relevant pathways in amphioxus and vertebrates (Pan et al., 2021; Curzon et al., 2023; Huang et al., 2023), we infer that wnt4a, foxl2 and cyp19a1 may play the major role in promoting ovarian differentiation of Asia arowana, while the expression and regulation of sox9, dmrt1, gsdf, amh and amhr2 possibly play an important role in development and differentiation of testis. Based on this inference, a potential regulatory network of Asian arowana sexual development was proposed (Figure 3B).

It is well known that TEs could regulate gene expression (Feschotte, 2008; Sundaram and Wysocka, 2020; Gebrie, 2023), and previous studies reported that the closer TE copy to the target gene, the higher correlation of their expressions (Dechaud et al., 2021). In mammalian genomes, TEs are an important source of various cis-regulatory sequences; 20% of the cis-regulatory elements (CREs) in the human genome may have been taken from TEs (Sundaram and Wysocka, 2020), and TEs often contribute to zebrafish cis-regulatory elements, tissue-specific expression and alternative promoters (Lee et al., 2022). Researchers have found that the fourth exon (Ex4) of the sex-determining gene dm-W in the African clawed frog (Xenopus laevis) originated from a non-coding fragment of the hAT-10 family of DNA transposons (Hayashi et al., 2022). An experimental evidence from sablefish (Anoplopoma fimbria) demonstrated that a TE insertion in the promoter region of gsdfY produced allelic diversification by bringing a cis-regulatory module, leading to transcriptional reprofiling and generating a new sex-determined gene for this species (Herpin et al., 2021). The regulatory elements of many genes contain TE sequences, which are involved in the regulation of gene expression. Some studies have reported that TEs can be expressed in multiple tissues (Brunet et al., 2018; Chang et al., 2022; Xu et al., 2023). However, because they align to multiple genomic loci, many RNA-seq reads derived from TEs are often discarded before data analysis. To resolve this problem, several computational tools have been developed (Jin et al., 2015; Jeong et al., 2017; Yang et al., 2019; Lanciano and Cristofari, 2020). We therefore observed a DNA transposon and an LTR transposon located at the promoter region of foxl2 expressed in the arowana ovary (Figure 4C), suggesting co-regulation of the cis-regulatory elements to enhance the neighboring sex-biased gene(s) and/or sex-biased TEs themselves on these two types of sequences, as well as the foxl2 may play a significant functional role in ovarian differentiation of arowana.

Dmrt1, an important member of the dmrt gene family, is essential for maintaining male-specific germ cells and testis differentiation (Matson and Zarkower, 2012; Cui et al., 2017). It not only participates in sex regulation in mammals, but also acts as a ubiquitous conserved sex regulation factor in other non-mammalian vertebrates including fishes. A recently published paper reported no orthologous gene of dmrt1 in amphioxus (Huang et al., 2023), while previous studies on the dmrt family of diverse animals found that dmrt1 and dmrtB1 are vertebrate-specific genes, and they may have arisen in vertebrates through WGD events (Mawaribuchi et al., 2019). Various studies on non-mammalian vertebrates have proved that dmrt1 or its paralogous genes control gonadal sex determination and differentiation with different mechanisms (Smith et al., 2009; Yoshimoto et al., 2010; Cui et al., 2017; Mawaribuchi et al., 2019). From the TF annotation of arowana transcriptomes, we identified multiple copies of dmrt, such as dmrt3 and dmrt2a. Through gene annotation and genome scanning of one elasmobranch and 18 bony fishes, we observed the conserved dmrt1-dmrt3-dmrt2 gene cluster among various vertebrates, and some species of the order Osteoglossus and Elopomorpha also have another dmrt3-dmrt2 gene cluster (Figure 5C). Previous studies have shown that dmrt1 may have emerged after a WGD event in vertebrates, and based on our findings we speculate that the dmrt gene clusters in these two basal teleost species are likely to come from the teleost-specific WGD event (3R-WGD).

It seems that the expanded dmrt gene cluster remains in this clade for continuous function after the 3R-WGD. This confirms that the two groups are sister branch to each other, and this branch forms as sister groups to all the other teleost (Parey et al., 2023). The gene cluster identified in the Asian arowana was named as dmrt1-dmrt3l-dmrt2l, based on the names of related genes published in the NCBI, and another cluster is named as dmrt3-dmrt2a (Figure 5C). Interestingly, dmrt1-dmrt3l in the dmrt1-dmrt3l-dmrt2l gene cluster was upregulated in the arowana testis (log2FC = 2.50, log2FC = 1.06), while dmrt2l was undetectable in both testes and ovaries; in the dmrt3-dmrt2a gene cluster, the dmrt3 gene is slightly expressed in the testes but without expression in the ovary, while dmrt2a is highly expressed in the ovary (log2FC = −3.39) (Figure 3A). These expression trends are very similar to those in other teleost fishes, indicating that dmrt1 and dmrt2 located in different clusters may be critical players in the sexual differentiation of Asian arowana: one of the gene clusters may be associated with male development, while another is potentially associated with female development.

Based on gonadal transcriptome sequencing and transcriptomic comparisons, we characterized a set of differentially expressed genes and pathways for potential involvement in sex determination or differentiation in Asia arowana. This sexual dimorphism is valuable for ascertaining the differential expression patterns of sex-related genes and enriching the genetic resources of this ornamental fish species. Our transcriptome data on sex-related genes can also promote the exploration of molecular mechanisms of gonadal development and sex determination in Asian arowana, which provides valuable references for practical gender identification and one-to-one breeding programs so as to expand fish number for a rapid elevation of worldwide economic value.