FOXL2 and DMRT1L Are Yin and Yang Genes for Determining Timing of Sex Differentiation in the Bivalve Mollusk Patinopecten yessoensis

Sex determination and differentiation have long been a research hotspot in metazoans. However, little is known about when and how sex differentiation occurs in most mollusks. In this study, we conducted a combined morphological and molecular study on sex differentiation in the Yesso scallop Patinopecten yessoensis. Histological examination on gonads from 5- to 13-month-old juveniles revealed that the morphological sex differentiation occurred at 10 months of age. To determine the onset of molecular sex differentiation, molecular markers were screened for early identification of sex. The gonadal expression profiles of eight candidate genes for sex determination or differentiation showed that only two genes displayed sexually dimorphic expression, with FOXL2 being abundant in ovaries and DMRT1L in testes. In situ hybridization revealed that both of them were detected in germ cells and follicle cells. We therefore developed LOG10(DMRT1L/FOXL2) for scallop sex identification and confirmed its feasibility in differentiated individuals. By tracing its changes in 5- to 13-month-old juveniles, molecular sex differentiation time was determined: some scallops differentiate early in September when they are 7 months old, and some do late in December when they are 10 months old. Two kinds of coexpression patterns were found between FOXL2 and DMRT1L: expected antagonism after differentiation and unexpected coordination before differentiation. Our results revealed that scallop sex differentiation co-occurs with the formation of follicles, and molecular sex differentiation is established prior to morphological sex differentiation. Our study will assist in a better understanding of the molecular mechanism underlying bivalve sex differentiation.


INTRODUCTION
Sexual reproduction is one of the most universal phenomena that widely exist in animals. As a focus of this research area, sexual development encompasses sex determination and differentiation. The former is defined as the process by which sex is established under the influence of genetic or environmental factors and the latter as the process that an undifferentiated gonad uses to transform into an ovary or a testis (Bull, 1983;Penman and Piferrer, 2008;Piferrer and Guiguen, 2008). Therefore, sexual development is a complex network that is initiated by a sex-determining trigger mediating the expression of sex differentiation genes, which ultimately gives rise to the phenotypic differences between sexes (Heule et al., 2014).
Sex determination and differentiation have received much research attention in various species. For example, in mammals, genotypic sex is defined by the presence or absence of sex-specific chromosome Y that carries the dominant male determinant SRY (Swain and Lovell-Badge, 1999;Kocer et al., 2009). In chickens that possess a ZZ/ZW sex chromosome system, DMRT1 is required for testis determination (Zarkower, 2001;Smith et al., 2009). In Caenorhabditis elegans and Drosophila melanogaster, a cascade of sex switch genes is controlled by the ratio of X chromosomes to sets of autosomes (the X:A ratio; Zarkower, 2001;Goodwin and Ellis, 2002). It seems sex is determined by sex chromosomes and controlled by some master switch genes in many organisms.
Mollusca represents the second largest phylum of invertebrates after Arthropoda. Bivalves are a large group of mollusks that exhibit different reproductive strategies: some are gonochoric, some are hermaphroditic, and some are capable of sex changes. Although much research has been done on many different bivalves, there is no clear evidence for the existence of sex chromosomes in these animals. However, some progress has been made regarding the identification of sex determination or differentiation genes. Several key sex-related genes in model species have been characterized in bivalves, including FOXL2 (Naimi et al., 2009a;Liu et al., 2012), WNT4 (Li et al., 2013;Yang et al., 2015), FST (follistatin; Ni et al., 2012), β-catenin (Li et al., 2014b;Santerre et al., 2014), DMRT (Naimi et al., 2009b;Feng et al., 2010;Yu et al., 2011;Shi et al., 2014), DAX1 (Li et al., 2014a), and SOXE (Santerre et al., 2014). Many of them display sexually dimorphic expression patterns in the gonads and are regarded as participants in the sex differentiation cascade. Meanwhile, gonadal transcriptome analyses identified some key candidate genes for sex determination or differentiation, such as FOXL2, WNT4, β-catenin, DMRT, DAX1, SOXE, and SOXH (Teaniniuraitemoana et al., 2014;Zhang et al., 2014;Tong et al., 2015;Li et al., 2016;Patnaik et al., 2016). These studies suggest that sex determination and differentiation genes may be deeply conserved in animals.
Yesso scallop Patinopecten yessoensis is a commercially important species widely cultured in China and Japan. It is predominantly gonochoric, with scarce hermaphroditism. Due to its commercial importance, much work has been performed to obtain an understanding of its reproductive process (Osada et al., 2004;Tanabe et al., 2010;Nagasawa et al., 2015). Recently, our group completed whole genome sequencing and gonadal transcriptome analysis of P. yessoensis (Li et al., 2016;Wang et al., 2017), which provide valuable resources for unraveling the molecular mechanisms underlying scallop sex differentiation. In the present study, morphological and molecular sex differentiation was examined in the Yesso scallop. Eight candidate genes for sex determination or differentiation were chosen to establish markers for sex identification and the onset of molecular sex differentiation was determined; we also compared the molecular changes in the undifferentiated and differentiated gonads. This study will pave the way for a better understanding of the regulatory network in bivalve sex differentiation.

Sample Collection
To obtain adult scallops with known sexes, 15-month-old mature female and male individuals were selected in May 2015, cultured in separate cages, and transported to the laboratory every month from August 2015 to March 2016. Meanwhile, juvenile scallops 5-13 months of age were obtained from July 2015 to March 2016. All of these scallops were collected from the Dalian Zhangzidao Fishery Group Corporation (Liaoning Province, China) and acclimated in filtered and aerated seawater for 1 week at the temperature at which they were collected. For each month, gonads of 50 adults and/or 50 juveniles were dissected: the majority were immediately frozen in liquid nitrogen and stored at −80 • C for RNA extraction, and the remainder were fixed in 4% paraformaldehyde overnight, dehydrated with serial methanol (25, 50, 75, and 100%) diluted in 0.01 M phosphate-buffered saline and then stored at −20 • C for paraffin sectioning and in situ hybridization. Our experiments were conducted according to the guidelines and regulations established by the Ocean University of China and the local government.

Histology
The samples were transferred to ethanol, cleared in xylene, embedded in paraffin wax, and cut into 5-µm-thick sections on a rotary microtome (Leica, Wetzlar, Germany). Serial sections were tiled on glass slides, deparaffinized with xylene, hydrated with graded ethanol to water, and stained with hematoxylin. After that, the glass slides were counterstained with eosin, dehydrated with ethanol, cleared with xylene, mounted with neutral balsam, and covered with coverslips. Finally, the sections were observed under a Nikon's Eclipse E600 research microscope.

Candidate Genes Identification
Four ovary-related genes (FOXL2, WNT4, FST, β-catenin) and four testis-related genes (DMRT1L, DAX1, SOXE, and SOXH) were selected for analysis. To identify these genes, protein sequences from other organisms were collected from NCBI and used as queries against the genome (Wang et al., 2017) and gonadal transcriptomes (Li et al., 2016) by TBLASTN with an e-value threshold of 1e −5 . The resultant scallop sequences were further confirmed by BLASTX against the NCBI protein sequence database.

Phylogenetic Analyses
Full-length protein sequences encoding FOXL2/3 and DMRT homologs were downloaded from NCBI. Multiple sequence alignments (Supplementary Material) were conducted using ClustalW2 program with default parameters (Larkin et al., 2007). The neighbor-joining phylogenetic trees were constructed using MEGA 6 (Tamura et al., 2013), with a Poisson model and assuming uniform rates among sites. Bootstrapping with 1000 replications was conducted to evaluate the phylogenetic tree.

RNA Extraction
Total RNA was isolated using the conventional guanidinium isothiocyanate method and digested with DNase I (TaKaRa, Shiga, Japan) to remove potential DNA contamination. RNA concentration and purity were determined by Nanovue Plus spectrophotometer (GE Healthcare, Piscataway, NJ, United States), and RNA integrity was verified by agarose gel electrophoresis. Only RNA samples with clear bands corresponding to 18S and 28S rRNA on the gel, an OD260/OD280 ratio between 1.8-2.0, and an OD260/OD230 ratio higher than 2.0 were used for subsequent experiments.

Reverse Transcription Quantitative-PCR (RT-qPCR)
First-strand cDNA was synthesized from 2 µg total RNA using oligo(dT) 18 and MMLV reverse transcriptase (TaKaRa, Shiga, Japan) in a volume of 20 µl. The reaction was performed at 42 • C for 90 min and terminated by heating at 70 • C for 10 min. Finally, the cDNA products were diluted to 10 ng/µl and stored at −20 • C. Gene-specific primers were designed using Primer Premier 5.0 and listed in Table 1. Amplification efficiency of each primer pair was calculated based on the standard curve generated from a twofold dilution series spanning six orders of magnitude. Quantitative-PCR was conducted using Light Cycler 480 SYBR Green I Master on a Light Cycler 480 Real-time PCR System (Roche Diagnostics, Mannheim, Germany) with the following program: 94 • C for 10 min, followed by 40 cycles of 94 • C for 15 s and 60 • C for 1 min. For each month, six to eight samples were assayed and all reactions were conducted in triplicate. To ensure that the RT-qPCR ran properly, negative controls, including NTCs (no-template controls) and no-reverse transcription controls, and positive controls were included in each run. Melting curve analysis was performed to verify that each primer set amplified a single product (Supplementary Figure S1). EF1A (elongation factor 1-alpha), which was stably expressed throughout the entire experiment, was used as an endogenous control for the normalization of gene expression (Santerre et al., 2013;Teaniniuraitemoana et al., 2015;Li et al., 2016). The relative expression level of each gene was calculated using the 2 − Ct method. Statistical analysis was performed by Welch's t-tests. Pearson's correlation coefficient was calculated to explore the relationships between the two genes. P-values lower than 0.05 were considered statistically significant.

In situ Hybridization
To prevent cross-detection of other FOX and DMRT genes, cDNA fragments avoiding conserved DNA binding domains were amplified with specific primers ( Table 2) containing a 5 T7 promoter sequence (5 -TAATACGACTCACTATAGGG-3 ). Purified PCR products were used as templates for in vitro transcription. Digoxigenin-labeled sense and anti-sense probes were generated using the DIG RNA Labeling Mix (Roche, Mannheim, Germany) and T7 RNA polymerase (Thermo, Waltham, MA, United States). Sections of the gonadal tissues were serially rehydrated in PBST (phosphate-buffered saline plus 0.1% Tween-20) and digested with 2 µg/ml proteinase K at 37 • C for 15 min. After pre-hybridization at 60 • C for 4 h, hybridization was performed with 1 µg/ml denatured RNA probe in hybridization buffer (50% formamide, 5× SSC, 100 µg/ml yeast tRNA, 1.5% blocking reagent, 0.1% Tween-20) at 60 • C for 16 h. Then, the probes were washed away, and antibody incubation was performed in a fresh solution of anti-digoxigenin-AP Fab fragments (Roche, Mannheim, Germany) coupled with blocking buffer (diluted 1:2000) at 4 • C for 16 h. After extensive  Primers F and R were used to amplify the cDNA fragment, and the product was used for the second-round PCR with the primers F and R-T7 to generate an antisense probe and with primers R and F-T7 to generate a sense probe.

Histological Analysis of the Juvenile Gonads
Gonads were examined histologically in juveniles aged 5 (shell height ∼10 mm) to 13 months old (shell height ∼60 mm). Figures 1A,B show the morphology of 5-and 6month-old gonads, respectively. As seen, the majority of the gonads are intestine, surrounded by gonadal tissues. In the 7-month-old gonad, some small-sized follicles have formed ( Figure 1C). In the 8-month-old gonad, the follicles grow bigger but are basically empty, containing some follicle cells and sexually indistinguishable gonia ( Figure 1D). The 10-month-old gonads are already morphologically differentiated, with oocytes ( Figure 1E) or spermatocytes ( Figure 1F) scattering in the follicles. Oogonia, oocytes, and a few mature oocytes coexist in the 11-month-old ovary (Figure 1G), and several layers of spermatogonia and spermatocytes are found in the 11-monthold testis (Figure 1H), indicating that the gonads are at the growing stage. As shown in Figures 1I,J, the 13-month-old gonads have reached full maturity. In the ovary, the follicles are filled with mature oocytes with a polygonal shape due to packing ( Figure 1I). In the follicles of the testis, diverse germ cells can be detected, including spermatogonia, spermatocytes, spermatids, and spermatozoa ( Figure 1J).

Selection of Candidate Genes for Sex Analysis
To determine the onset of molecular sex differentiation, molecular markers were screened for sex identification. All eight candidate genes (FOXL2, WNT4, FST, β-catenin, DMRT1L, DAX1, SOXE, and SOXH) were identified in the scallop genome. Their expression patterns in ovaries and testes were examined by RT-qPCR. In order to obtain genes that show sexually dimorphic expression in somatic and/or gonia cells and eliminate the effects of gametogenesis, we only used gonads at the resting stage for gene expression analysis (Yu et al., 2017). As shown in Figure 2A, only FOXL2 and DMRT1L displayed significant (P < 0.01) sexually dimorphic expression, with FOXL2 being abundant in ovaries, and DMRT1L in testes. Therefore, LOG 10 (DMRT1L/FOXL2) is a potential marker for scallop sex identification, which led us to further investigate its values in the other three reproductive stages (proliferative, growing, and maturation stages). The results also confirmed that LOG 10 (DMRT1L/FOXL2) could be used to determine the sex of differentiated scallops throughout the reproductive cycle: the  . " * * " indicates differences that are statistically significant (P < 0.01). (B) The LOG 10 (DMRT1L/FOXL2) in ovaries (orange) and testes (blue) at the resting stage, proliferative stage, growing stage, and maturation stage. For each stage, six gonads, including three ovaries and three testes, were assayed. The horizontal lines indicate the threshold for sex differentiation, i.e., orange for ovary and blue for testis, respectively.
value was always lower than 0 for ovaries and higher than 2 for testes ( Figure 2B).

Phylogenetic Analyses of FOXL2 and DMRT1L
To determine the identity of the two sexually dimorphic expressed genes, phylogenetic analysis was conducted. As shown in Figure 3A, Yesso scallop FOXL2 clustered with FOXL2 from other bivalves and gastropods, and then clustered with vertebrate FOXL2 and FOXL3. In Figure 3B, the selected DMRT proteins were clustered into several major groups. The Yesso scallop DMRT, together with DMRT from Mimachlamys nobilis and Haliotis asinina forms a new group, namely, DMRT1L.
Localization of FOXL2 and DMRT1L mRNA Spatial expression of FOXL2 and DMRT1L was detected in the ovary and testis at the proliferative stage after sex differentiation. For both genes, the signals were present in germ cells and follicle cells. Specifically, the anti-sense probe of FOXL2 was detected in the cytoplasm of oogonia, oocytes, and follicle cells ( Figure 4A) and faint in testis ( Figure 4B). DMRT1L transcripts were located clearly in all germ cells of the ovary ( Figure 4C) and testis ( Figure 4D), but the signal intensity was higher in the testis than in the ovary. No signal was detected in the ovary or testis with the sense probes (Figures 4E-H).

Determining the Timing of Molecular Sex Differentiation
The feasibility of LOG 10 (DMRT1L/FOXL2) for accurate sex identification enables us to determine the timing of gonadal sex differentiation in the Yesso scallop. Here, we traced the dynamic changes of LOG 10 (DMRT1L/FOXL2) in juveniles aged 5-13 months old. According to Figure 5A, LOG 10 (DMRT1L/FOXL2) values of young individuals (5 and 6 months old) are between 0 and 1, suggesting that these scallops are still undifferentiated.
In 7-to 10-month-old scallops, we found a decreasing number of LOG 10 (DMRT1L/FOXL2) falling within 0-1. This indicates that sex differentiation does not occur simultaneously in all individuals, with some scallops differentiating early when they are 7 months old and some not differentiating until they are 10 months old. However, after reaching 11 months of age, all investigated gonads had differentiated into ovaries or testes. The histogram in Figure 5B shows the corresponding differentiation rate for each month based on the LOG 10 (DMRT1L/FOXL2) values. It reveals an increase in the differentiation rate from 0 in 5-and 6-month-old individuals to 100% in 11-to 13-month-old scallops.
Coexpression of FOXL2 and DMRT1L in the Juvenile Gonads Figure 6A shows the temporal expression profiles of both FOXL2 and DMRT1L in juveniles of 5-to 13-month-olds. According to the results, FOXL2 and DMRT1L are relatively low in 5-to 10-month-old scallops, but expression of FOXL2 in ovaries and DMRT1L in testes was dramatically increased in 11-to 13-monthold scallops. During the eight investigated months, there existed two major coexpression patterns between FOXL2 and DMRT1L: (i) a significantly positive correlation for the two undifferentiated stages (r = 0.955-0.980, P < 0.01), with the regression equation of y = 1.018x + 0.370 ( Figure 6B); (ii) significantly negative correlation after the completion of sex differentiation in the last three time points (r = −0.791 to −0.957, P < 0.05), with the regression equation of y = −0.684x + 1.584 ( Figure 6C). Between the undifferentiated and fully differentiated months are three partially differentiated ones, in which correlation coefficients displayed a transition from significantly positive (0.860) to negative (−0.587).

DISCUSSION
By examining the gonadal expression of eight candidate sex determination or differentiation genes in the resting stage, we found that only two genes (FOXL2 and DMRT1L) displayed sexually dimorphic expression. Both genes are transcription factors, with the former being a forkhead transcription factor and the latter encoding a DM domain-containing protein. Specific sequence information, including conserved domain alignments, has been described previously (Li et al., 2016). According to the transcriptomes of nine adult tissues/organs (Wang et al., 2017), FOXL2 and DMRT1L seem to be specifically expressed in scallop gonads, and no evidence of splice isoforms was found for DMRT1L despite that splicing of DMRT is common in many animals (Nothiger et al., 1987;Kopp, 2012). Interestingly, sexually dimorphic expression of FOXL2 (Liu et al., 2012;Teaniniuraitemoana et al., 2014;Zhang et al., 2014) and DMRT1L (Shi et al., 2014) occurs not only in the Yesso scallop but also in other bivalves. Therefore, FOXL2 and DMRT1L, like Yin and Yang, may be crucial genes for female and male gonadal differentiation in bivalves. Similarly, in mammals, FOXL2 is an ovarian marker that suppresses genes involved in testicular differentiation from the early embryonic gonad throughout adult life (Cocquet et al., 2002;Ottolenghi et al., 2007;Uhlenhaut et al., 2009;Liu et al., 2014). Although DMRT genes have diverse roles in development and physiology (Huang et al., 2005;Hong et al., 2007;Reitzel et al., 2016), DMRT1 is a deeply conserved gene that is involved in testis differentiation and development in various animal taxa (Naimi et al., 2009b;Kopp, 2012;Matson and Zarkower, 2012;Picard et al., 2015). In mammals, a battle between FOXL2 and DMRT1 for sex control has been proposed (Herpin and Schartl, 2011;Matson et al., 2011). We assume that this battle may be extended to invertebrates, but with different DMRT members. Spatial expression of FOXL2 and DMRT1L in the ovary and testis reveals that both genes are distributed in germ cells and follicle cells, consistent with previous research in other marine bivalves (Naimi et al., 2009a;Liu et al., 2012). The signal intensity also agrees with the RT-qPCR data. The expression patterns of FOXL2 and DMRT1L make LOG 10 (DMRT1L/FOXL2) an alternative method to histological examination for scallop sex identification. By tracking LOG 10 (DMRT1L/FOXL2) at four gonadal developmental stages, we found the value was always lower than 0 for ovaries and higher than 2 for testes, which indicates that in the Yesso scallop, (i) there is no sex reversal once the gonads differentiate, (ii) FOXL2 and DMRT1L may play important roles in both sex maintenance and gametogenesis, and (iii) LOG 10 (DMRT1L/FOXL2) is an effective solution for sex identification during the resting stage, which is usually indistinguishable by morphology-or histologybased approaches. Since FOXL2 and DMRT1L also show sexually dimorphic expression patterns in other bivalves (Liu et al., 2012;Shi et al., 2014;Teaniniuraitemoana et al., 2014;Zhang et al., 2014), using LOG 10 (DMRT1L/FOXL2) for sex identification may be applicable in these organisms. For oysters that are capable of sex reversal (Teaniniuraitemoana et al., 2014;Zhang et al., 2014), sex-reversing individuals may be easily screened out by examining the LOG 10 (DMRT1L/FOXL2) values that fall between    Figure 5A. The numbers above the bars are the correlation coefficients between LOG 10 (FOXL2) and LOG 10 (DMRT1L) in the particular months. " * " indicates P < 0.05; " * * " indicates P < 0.01. (B) The positive correlation between LOG 10 (FOXL2) and LOG 10 (DMRT1L) for 5-to 6-month-old undifferentiated gonads. (C) The negative correlation between LOG 10 (FOXL2) and LOG 10 (DMRT1L) for 11-to 13-month-old differentiated gonads.
the ovary and testis thresholds, which could contribute to a better understanding of the molecular mechanisms underlying sex reversal.
Based on the histological analysis, we found that the developmental process of juvenile gonads can be classified into four stages, namely, undifferentiated, differentiated, growing, and maturation stages, which are consistent with the documentation by Shumway and Parsons (2016). Combined with the dynamic changes of LOG 10 (DMRT1L/FOXL2) in juvenile scallops, we found that sex differentiation occurs earlier at the molecular level than at the histological level, which is similar to the findings in other organisms (Siegfried, 2010;Haugen et al., 2012;Ayers et al., 2013;Robledo et al., 2015;Tao et al., 2018). In addition, our study shows that although the onset of molecular sex differentiation varies between individuals, it always occurs from 7 to 10 months of age, and almost all of the scallops investigated completed differentiation before 11 months of age. This falls within the time frame of gonadal differentiation documented by Shumway and Parsons (2016) (between 4 months old and 1 year old). It suggests that LOG 10 (DMRT1L/FOXL2) allows early identification of sex, and is a convenient and accurate way to assess the molecular sex differentiation period. This opens a gate for unraveling the regulatory network of sex differentiation in the Yesso scallop.
The coexpression patterns of FOXL2 and DMRT1L are different before and after gonadal differentiation. Before differentiation, FOXL2 and DMRT1L are relatively low but display unexpected coordinated expression patterns. Combined with the histological results showing that the follicles have not yet formed at this stage, we assume these two genes may coexpress in GSCs that are dispersed in the connective tissues. Sex differentiation co-occurs with the formation of follicles, in which relative expression of FOXL2 vs DMRT1L changes, possibly leading to the differentiation of undifferentiated germ cells into spermatogonia or oogonia. Therefore, we speculate that in Yesso scallops, the undifferentiated gonads directly differentiate into either ovaries or testes. This challenges a previous opinion regarding the existence of sex reversal from male to female through an undifferentiated stage (Maru, 1978). In particular, some researchers claimed that Yesso scallops grow as males at less than one year of age; sex is thereafter reversed to female through an undifferentiated or ambiguous status in half of the scallops, resulting in a sex ratio close to 1:1; after sex reversal, the sex may become stable throughout successive periods of life (Mori et al., 1977;Maru, 1978;Otani et al., 2017). Considering that these studies are generally based on macroscopic observation of gonadal color (testis with white color, ovary with orange or pink color), which could be misleading in some cases, we assume that findings of sex reversal in Yesso scallops may be an artifact. However, since environmental factors have been demonstrated to be involved in the sex determination of many bivalves, we cannot rule out the possibility that the Yesso scallop does undergo sex reversal in Japan.

AUTHOR CONTRIBUTIONS
LinZ and SW conceived and designed the experiments. RL, WL, and YZ performed the experiments. YL, MZ, LiaZ, and XH collected the samples. RL and LinZ analyzed the data. RL, LinZ, SW, and ZB wrote the paper. All authors read and approved the final manuscript.