The sample collection and animal experiments were conducted in accordance with the regulations of the Guide for Care and Use of Laboratory Animals and were approved by the Animal Care and Use Committee of Jimei University.

To obtain samples for RNA-seq, artificially bred yellow drum of 26 days post-hatch (dph) purchased from a breeding farm in Zhangzhou, Fujian province, China. Fries were reared in two tanks (80 cm × 60 cm × 60 cm) with temperature at 28°C and salinity at 25 ~ 30, dissolved oxygen at 6 ppm. Fries were fed with compound feed three times a day. The sample collection at five critical stages (30 dph, 36 dph, 46 dph, 55 dph, and 61 dph) of yellow drum gonadal development according to Lou et al. (Lou et al., 2016), approximately 800 individuals were collected at each stage. Fries were anesthetized with tricaine methanesulfonate (MS-222, Shanghai Acmec Biochemical Co., China) at 300 mg/L before sampling, and the body length was recorded (Table 1). The head and tail of each fry was pruned, and the rest of the body was cut open from the abdomen and stored in RNA Later solution. Subsequently, gonads were dissected and collected from the body under a microscope for RNA extraction and the remaining part was stored at 75% ethanol for DNA extraction. The sex of each sampled fry was determined through PCR amplification with the sex specific markers (Sun et al., 2018). To obtain a sufficient quantity of RNA for sequencing, we generated three pools of gonads isolated from the fries for each sex at each stage. Total RNA was isolated using TRIZOL Reagent (Invitrogen, USA) according to the manufacturer’s instructions. The purity and integrity of the RNA were verified with agarose gel electrophoresis and Bioanalyzer 2100 (Agilent Technologies, USA). RNA integrity number (RIN) ≥ 8.0 were subjected to RNA enrichment with magnetic beads with Oligo (dT), and RNA was disrupted into short fragments by adding fragmentation buffers. These short fragments were used as templates and random hexamers as primers to synthesize the first strand of cDNA, followed by the synthesis of second strand cDNA. The AMPure XP beads were employed to purify double-strand cDNA, EB buffer was used for end-repair and A-tailing. The purified fragments were ligated with Illumina sequencing adapters, and paired-end libraries were sequenced on an Illumina NovaSeq 6000 system.

The clean reads were obtained from raw data by filtering out adaptor contaminated reads and low-quality reads. The quality of the sequenced reads was assessed by FastQC v0.11.9 (Andrews, 2010). Clean reads were aligned to the male reference genome using STAR v2.5.3a (Dobin et al., 2013). After alignment, the transcripts were assembled using StringTie v1.3.3b (Pertea et al., 2015), and the number of reads that mapped to the genes and exons was calculated using FeatureCounts v2.0.1. The analysis of differentially expressed genes was performed with DESeq2, and the related R packages were employed for the PCA analysis. Significantly differentially expressed genes (DEG) were defined as adjust P-value < 0.05 and log₁₀2(Foldchange) > 1. Gene Ontology (GO) terms analysis of DEGs and Kyoto Encyclopaedia of Genes and Genomes (KEGG) pathway enrichment were performed using clusterProfiler R packages.

The interaction network of differential genes between females and males at 36 dph was analyzed using the STRING data [STRING: functional protein association networks (string-db.org)]. The CytoScape software was used to display the result.

The protein sequences of 11 species including Homo sapiens (human), Gallus gallus (chicken), Lctalurus punctatus (spotted gar), Danio retio (zebrafish), Oryzia latipes (medaka), Takifugu rubripes (fugu), Tetraodon nigrovirids (spotted green pufferfish), Oreochromis niloticus (tilapia), Gasterosteus aculeatus (three spined stickleback), and Larimichthys crocea (large yellow croaker) were obtained from Ensembl website (https://asia.ensembl.org/index.html). The phylogenetic tree was constructed with MEGA 5.0 using UPGMA method. The adjacent genes with dmrt1 of the 11 species are referred to the Ensembl Genome Browser (http://www.ensembl.org).

The Kank1 gene reference sequences were extracted from transcriptome data and the amino acid sequences were translated through DNAMAN. The nuclear localization signal (NLS) was predicted using online tools NLStradamus (utoronto.ca) and the nuclear export signal (NES) was predicted using online tools LocNES NES prediction tool provided by Chook Lab (swmed.edu). The Kank1 gene expression was assessed using FPKM at different development stages in yellow drum and significant differences between the males and females were determined by t-test (p < 0.05).

Kank1 expression vector (pCS2-dmrt1-mCherry) was constructed and the coding sequence was obtained through PCR cloning from cDNA library of yellow drum. The pCS2- kank1-mCherry vector was double digested with BmaHI and ClaI, then pCS2-dmrt1-mCherry was constructed with dmrt1 coding sequence cloned into the digested site using ClonExpress^® II One Step Cloning Kit (Vazyme, Nanjing, China). In order to distinguish with mCherry fluorescence, we amplified the EGFP coding sequences using the pEGFP-N1 vector as a template and replaced the mCherry coding sequences of the vector pCS2-dmrt1-mCherry using the ClaI and XhoI double digestion to construct the vector pCS2-dmrt1-EGFP. The dmrt1 promoter sequences were amplified from the genomic DNA of the yellow drum, the product was double digested with KpnI and HindIII, then cloned into the same site of pGL3-basic vector (Promega, USA) to construct the vector pGL3-basic-dmrt1_promoter. All the primers used in PCR cloning warelisted in Table 2.

CHO (a cell line from Cricetulus griseus ovary) and HEK293 were cultured in 1640 and high glucose DMEM medium (Dulbecco’s modified eagle medium, Gibco), respectively, supplemented with 10% fetal bovine serum and 1% penicillin 100 units/ml with 100 ug/ml streptomycin, plated in 25 cm² flask, and placed in the carbon dioxide incubator. The detailed operation steps were followed according to the reference manual.

Lipo8000™ Transfection Reagent was used in vector transfection (Beyotime Biotechnology, China). Cells were plated in 12-well plates at 18 ~ 24 hours before transfection. While the cells were completely adherent and stuck to the wall, the fresh medium was replaced according to protocol. The cells were transfected with 50 μL Opti-MEM Medium, 1 μg DNA and 1.6 μL Lipo8000™. The ratio of the pGL3-basic promoter vector and the pRT-TK vector (an internal control) was about 40: 1. In the co-transfection assay, the ratio of expression vector to promoter vector was 1: 48, 2: 48, 4: 48, 8: 48, respectively. Luciferase activity was measured at 24 hours after transfection using the dual-luciferase reporter assay system (Beyotime Biotechnology, China) and LB 9508 Lumat³ (Berthold technologies, Germany) following the manufacturer’s protocol. Assays were performed in triplicate and expressed as mean ± standard error.

Firstly, the medium in the 12-well cell plates was discarded and washed once with 1 ml PBS. Cells were then fixed with 1 mL 4% paraformaldehyde for 10 minutes, and were permeabilized with 0.2% Tween 100 for 10 minutes and washed once again with 1 mL PBS. Nuclear staining was performed with 1 mL 0.2% DAPI for 15 minutes, and observation and photography were done with a laser confocal scanning microscope (TCSSP8, Leica).

RNA libraries generating from the gonad of the various early stages were sequenced on an Illumina NovaSeq 6000 system. The total raw data and clean data were 318.7 Gb and 307.7 Gb, respectively. The average reads map rate was 87.9% when aligned the clean reads to the assembly of yellow drum with STAR. Transcriptome assembly obtained 23,997 transcripts from 18,596 genes (Table 3).

We performed PCA analysis with all the samples (Figure 1). Clustering discovered that the samples were similarly grouped into two sections by the sex (dimension 1, Figure 1A) and the stage (dimension 2, Figure 1B). This result suggests that the samples before 36 dph formed a group and the samples after 36 dph formed another group, indicated 36 dph (body length = 24.4 ± 2.2 mm) was a critical stage. Subsequently, we analyzed the different expression genes between females and males at each stage (Figure 2A) and counted the genes of up- or down-regulated (Figure 2B), and found that the number of genes with up- and down-regulated at 36 dph were more than at other stages. This could be due to the sex biased genes expressed in two kinds of gonad before 36 dph, and switched into sex-specific expression patterns once the sex-determining key gene was expressed. A previous study showed that the critical period of sex differentiation in the yellow drum occurred at the stage of 36 dph (total length = 21.2 ± 3.9 mm) (Lou et al., 2016). Our results imply the 36 dph is a sex determination critical stage, which is inconsistent with the previous result.

To investigate the function of these differential genes, we performed GO and KEGG analyses between females and males at 36 dph (Figure 3). In females, the upregulated genes mainly participated in the biological processes associated with development process, reproductive process, gamete generation, germ cell development, and cellular process involved in reproduction in multicellular organisms (Figure 3A). Of note, the cell cycle and Wnt signaling pathway were primarily enriched (Figure 3C). In males, upregulated genes were mainly involved in transmembrane receptor protein tyrosine kinase signaling pathway, cell migration, cell motility, and tube morphogenesis (Figure 3B). The KEGG pathway analysis suggested that Focal adhesion and TGF-beta signaling pathway played important roles in males at 36 dph (Figure 3D). According to the general process of sexual development, we speculated that the primordial germ cells had migrated to the place of primitive gonad and began to proliferate and develop in females. However, the gonad development was relatively late in males. The primordial germ cells may be undergoing the cell migration process. In addition, we collected genes that have been reported to express in germ cells, somatic cells (Matsuda et al., 2002; Richardson and Lehmann, 2010). Comparing the genes expressed in germ cells, somatic cells, and the genes involved in the biological processes of cell cycle and cell migration at different stages (Figure 4), we found that the genes participating in the cell migration process expressed higher in males than in females. However, the genes expressed in germ cells involved in the cell cycle expressed higher in females. The gene expression trend was consistent with the results of GO and KEGG enrichment.

The differential genes pathway analysis provided a clue that mabe the different behavior of the primordial germ cells caused the differences in gonad development between males and females. Given dmrt1 plays a pivotal role in the sex determination of yellow drum, we wanted to understand the link between master sex determining gene Dmrt1 and the development of primordial germ cells. We analyzed the protein-protein interaction networks with the differential genes associated with the behavior of primordial germ cells from the pathway analysis (Figure 5). The result implies that the proteins interaction mainly involved in the cell motility process and in the cell cycle in males and females, respectively. The network of proteins interaction was complex. We paid close attention to which proteins may interact with Dmrt1 in yellow drum, including Sox8, Sox9, Fgf6, Bmpr2, Kank1, Foxl2, Wnt4a, and Orc1. However, it seemed that only Kank1 could act upon Dmrt1 directly among the various genes. This result indicated that the Kank1 might take part in the germ cells migration process and could interact with Dmrt1 in yellow drum. Previous studies reported that Kank1 negatively regulated actin polymerization and cell motility through inhibiting the binding of Rac1 and insulin receptor substrate p53 (IRSp53) (Mhawech, 2005; Kakinuma et al., 2008; Kakinuma et al., 2009; Roy et al., 2009). Therefore, we suspected that Kank1 might inhibit the germ cells migration at the 36 dph in male individuals of yellow drum. However, the relationship between Dmrt1 and Kank1 is unclear, and has not been reported in previous related studies.

We found that the kank1 gene os located in the sex chromosome near dmrt1. Previous studies have shown that sex chromosomes usually evolved from autosomes, and once the sex determining gene fixed in the chromosome, the sexually antagonistic allele genes would be located in the surroundings (Charlesworth et al., 2005). It seemed the kank1 gene may be the sexually antagonistic allele gene. We then constructed the phylogenetic trees of dmrt1, and analyzed the collinearity between kank1 and dmrt1 in different species (Figure 6). The result suggests that the collinearity between kank1 and dmrt1 was conserved from fishes to humans. This result suggested that the kank1 gene might act as an important role in the sex determination in yellow drum, and it should be very meaningful to study the gene function of kank1.

The ORF of kank1 consisted of 3993 bp encoding 1330 amino acid residues in yellow drum (Figure 7A). A KN motif and an Ank repeat were predicted in the C-terminal and N-terminal, respectively. The nuclear localization signal (NLS) and the nuclear export signal (NES) were predicted (Figure 7A), implying Kank1 may be a shuttling protein. The kank1 gene expressed significantly higher in males than females from 46 dph to 55 dph, which was the stage that the cytological characteristics of male gonad differentiation were gradually observed (Figure 7B). This result suggests that the kank1 gene might act as an indispensable role in the male gonad differentiation.

A previous study showed Dmrt1 was different in the X and Y chromosome, respectively named as Dmrt1X and Dmrt1Y. Dmrt1, as a transcription factor, can bind to its own promoter and regulate the expression of itself (submitted, Sun et al.). We constructed the expression vector of Dmrt1 and Kank1, and co-transfected them into the CHO cells. The result suggests that Dmrt1 expressed in cell nucleus and Kank1 expressed in karyoplasm (Figure 8), indicating Kank1 was really a shuttling protein. With the dual-luciferase reporter assay, we detected the transcription level of dmrt1_Luc gene was higher when transfecting the pCS2-dmrt1-EGFP and pCS2-kank1-mCherry with pGL3-basic-dmrt1_promoter simultaneously than only transfecting the pCS2-dmrt1-EGFP with pGL3-basic-dmrt1_promoter in CHO cells and 293T cells (Figures 9A, B). This result demonstrates that Kank1 together with Dmrt1 enhanced the expression of dmrt1, possibly through its interaction with Dmrt1 protein to act on the promoter of dmrt1. Previous work has shown that Kank1 negatively regulated actin polymerization and cell motility (Kakinuma et al., 2009). We speculated that Kank1 might take part in the sex determination through inhibiting the primordial germ cells migration at the cell development and increased the expression of sex-determining gene - dmrt1- in yellow drum.

Gonad is constituted by somatic cells and germ cells. And sex determination includes somatic sex determination and germ cell sex determination (Nagahama et al., 2021). It has long been believed that gonadal somatic cells instruct germ cells to adopt their sexual fates, and somatic signals decide germ cells whether they reside in a male or female body. The sex-determining genes that have been reported mostly express in somatic cells and the mutation of these genes result in sex reversal, such as Sry (Sinclair et al., 1990), Dmy (Matsuda et al., 2002), and Amhy (Hattori et al., 2012). Sex reversal of genetic females (XX) could be induced by the transplantation of XY somatic cells in the medaka, Oryzias latipes (Shinomiya et al., 2002). Nishimura et al. found that foxl3 was a germ cell-intrinsic factor involved in sperm-egg fate decision in medaka. XX medaka with disrupted foxl3 developed functional sperm in the expanded germinal epithelium of a histologically functional ovary (Nishimura et al., 2015). A recent report showed that the sexual fate of tilapia germline was determined by the antagonistic interaction of dmrt1 and foxl3 (Dai et al., 2021). The conventional idea is challenged by these recent studies. These results also indicate that the crosstalk between germ cells and gonadal somatic cells is critical for sex differentiation of the gonads. To address the question, the relevant research in this study was carried out.

In this study, we confirmed the 36 dph was a sex determination critical stage in yellow drum through transcriptome analysis at the early stages. Subsequently, we analyzed the differential genes and the signaling pathways they participated in at 36 dph. We discovered that the signaling pathways in females and males were different at sex determination critical stage. In females, the cell cycle and Wnt signaling pathway were mainly enriched and the most biological process was about germ cell development and gamete generation. This result suggests the germ cells in the differentiation process in females. Given germ cells confront two major cell fate decisions as they move from an immature state into sexuality, one decision is entering into a germline-specific cell cycle (i.e., meiosis), and the other is the commitment to differentiate as spermatozoa or eggs (Kimble and Page, 2007), and the cell cycle entering into meiosis is considered as a hallmark of feminization of germ cells (Kikuchi et al., 2020). We suspected that the germ cells sex determination would progress toward to the females when it was in cell cycle and reproduction. Dmrt1 has been cloned and was demonstrated as the most possible sex determination gene in yellow drum (Sun et al., submitted). In medaka, the dmrt1 gene resulted in negative regulation of germ cells proliferation and mediated a mitotic arrest of germ cells (Herpin et al., 2007). In yellow drum, the dmrt1 gene almost exclusively expressed in males, implying it was a barrier for proliferation and feminization of germ cells, and an intrigue for the gonad masculinization. In males, the main enrichment signaling pathway included TGF-beta signaling and Focal adhesion, the main biological process involved cell migration, cell motility, and tube morphogenesis. So, another possibility was the dmrt1 gene may negatively regulate the germ cells migrated to the primordial gonadial ridge and caused two types of gonadal tissues to be observed in the development stage of 40 days post-hatch (dph), which were designated as tissue type I and type II (Lou et al., 2016). Unfortunately, we failed to investigate the numbers of germ cells at early stage using histological experiments in yellow drum. Future research would pay more attention to this question.

We found the Kank1 protein contained a nuclear localization signal (NLS) and a nuclear export signal (NES), and expressed in karyoplasm in CHO cells. These results were consistent with the Kank protein in humans (Wang et al., 2006), and suggest that Kank1 of yellow drum was also a shuttling protein. Nuclei-cytoplasmic transportation of proteins plays an important role in the regulation of many cellular processes, including signal transduction, cell cycle progression and cell proliferation (Gama-Carvalho and Carmo-Fonseca, 2001). Furthermore, studies have demonstrated that Kank could bind to β-catenin and promote the nuclear localization of β-catenin, which results in the activation of β-catenin-mediated transcription (Wang et al., 2006). The β-catenin is an important component in the Wnt signaling pathway which is involved in various differentiation events during embryonic development (Lustig and Behrens, 2003). Evidence in mice and humans implied WNT4 and R-SPONDIN-1 (RSPO1) promoted the ovarian fate and blocked testis development, and activation of β-catenin was sufficient to block the testis pathway (Maatouk et al., 2008). In protandrous black porgy (Acanthopagrus schlegeli), wnt4 has important roles in late ovarian growth and natural sex change, such as oogonia proliferation, structure of ovarian lamellae, and vitellogenic oocytes (Wu and Chang, 2009). A previous study reported that rspo1 expressed more strongly in ovaries than in testes in zebrafish (Zhang et al., 2011). In Nile tilapia (Oreochromis niloticus), β-catenin1 and β-catenin2 played roles in ovarian differentiation and maintenance of female pathway (Wu et al., 2016). However, the roles of the RSPO1/Wnt/β-catenin signaling pathway in the sex determination in fish is still not clear. In teleost, the sequential hermaphroditism is not unusual (Colombo and Grandidr, 1996; Helfman et al., 1997), protogynous fishes develop into females first and become males later during fry stage, such as zebrafish (Maack and Segner, 2003), lager yellow croaker (You et al., 2012), and yellow drum (Lou et al., 2016). In the female-determination process of yellow drum, it is not clear whether Kank1 could activate the β-catenin-mediated transcription and exert the same function as in mice and humans, given that kank1 was not differentially expressed between genetic females and males before 36 dph.

The kank1 gene was expressed significantly higher in males than females at sex-determination critical stages, indicating kank1 was more favorable in male sex determination. With the protein network interaction analysis, we uncovered that Kank1 took part in the germ cells migration process and interacted with Dmrt1. Kank1 might act as a shuttling protein as in humans and play difunctional roles at the early gonadal development stages, i.e., not only participated in the germ cells motility, but also increased the expression of dmrt1 in yellow drum. Although the histological location of Kank1 has not been confirmed, we speculate it should at least express in germ cells based on the role of Kank1 in the cell migration in yellow drum. Kank1 enhanced the expression of dmrt1 in yellow drum, and dmrt1 is the master sex-determining gene and triggered the testis development. These results indicate kank1 might act as an upstream factor of the dmrt1 gene. Of note, the sex battle and antagonism between female and male are from beginning to end. So, the balance of these factors leading to either male or female development is important. Once the balance is destroyed, the sex reversal will happen.