Identification of a goat intersexuality-associated novel variant through genome-wide resequencing and hi-c

Background

In as early as the nineteenth century, people regarded hornlessness as a beneficial and important economic trait and bred specialized hornless goat strains. However, during breeding, the proportion of intersex individuals in the hornless goat population gradually increased. This phenomenon was termed polled intersex syndrome (PIS) ( Eaton, 1945 ). Intersexuality, the phenomenon wherein certain dioecious organisms possess both sexes, has been widely observed in various livestock species ( Bosu and Basrur, 1984 ; Wang and Zhang, 1993 ), including goats ( Ramadan and El Hassan, 1988 ; Ramadan et al., 1991 ), within the last century. The proportion of intersex goats within the global population is 4–15% ( Zhan et al., 1994 ; Chen et al., 2010 ; Song et al., 2015 ). Reproductive system malformations in PIS goats lead to the loss of reproductive capacity and are thus some of the great challenges encountered in the development of the goat industry.

In 1996, the CHI1q41–q45 genomic regions were confirmed to be linked to hornlessness ( Vaiman et al., 1996 , 1999 ). Various molecular methods, such as chromosome walking technology and sequencing, have been used to refine the PIS locus to <100 kb ( Schibler et al., 2000 ). In 2001, the indicator of PIS in goats was mapped and resolved to an 11. 7 kb non-coding deletion in CHI1q43 that was located ~200 kb upstream of the FOXL2 gene ( Pailhoux et al., 2001 ). FOXL2 is an important sex determination gene ( Bagheri-Fam et al., 2017 ; Tao et al., 2018 ) with a key role in ovarian development ( Pannetier et al., 2016 ; Elzaiat et al., 2017 ). For example, a previous study on mouse models with ovarian FOXL2 gene deletion showed that FOXL2 +/− mice have a normal phenotype, FOXL2 +/+ mice have a similar phenotype as patients with human blepharophimosis syndrome, and FOXL2 −/− mice exhibit narrow eye slits and premature ovarian failure ( Baron et al., 2005 ). Furthermore, by using gene editing technology, Boulanger et al. (2014) verified that the loss of function and expression silencing of FOXL2 can cause female-to-male sex reversal in XX goats ( Boulanger et al., 2014 ).

Moreover, the development of PIS diagnostic molecular markers can effectively avoid the misdiagnosis caused by the phenotypic identification of obscure PIS cases. Although a series of PIS diagnostic methods based on PCR amplification has been reported ( Yang et al., 2012 ; Zhang et al., 2019 ), some studies on the diversification of PIS deletion structure have questioned the accuracy of these methods ( Li et al., 2011 ; Kijas et al., 2013 ). For example, some intersex Rangeland goats do not exhibit the known homozygous PIS deletion ( Kijas et al., 2013 ). Therefore, whether PIS deletion is specific for the diagnostics of intersexuality in goats remains controversial. Notably, the long-read whole-genome sequencing of two (one Saanen and one Valais Blacknecked black goats) genetically female (XX) intersex goats ( Simon et al., 2020 ) demonstrated that a highly complex structural variant involving a ~0. 48 Mb duplicated segment from ~21 Mb of chromosome 1 (CHI1) is inversely inserted into the known PIS deletion and that the length of the PIS deletion has also been shortened to 10. 159 kb from 11. 7 kb ( Pailhoux et al., 2001 ).

In this study, we, for the first time, identified the intersex-related genetic variation structure of the Chinese goat population via high-throughput sequencing technology and analyzed the chromosomal spatial structure of the PIS-related genetic structure through high-throughput chromosome conformation capture (Hi-C) technology to obtain an in-depth understanding of the molecular genetic mechanism of PIS. Our work could also provide a valuable reference for the future development of diagnostic tools with enhanced broad-spectrum recognition capabilities.

Methods

Genomic Library Construction and Sequencing

All the experimental conditions of this study were approved by the Committee on the Ethics of Animal Experiments of the Southwest University (No. [2007] 3) and the Animal Protection Law of China.

We collected venous blood samples from 55 goats comprising 35 intersex goats (26 XX Tangshan dairy goats and 9 XX Chinese southern native goats) and 20 XX non-intersex Tangshan dairy goats. A total of 2 mL venous blood was collected from each animal (Sampling from Tangshan dairy goat breeding farm, Tangshan, China). The wound was sterilized with 70% medical alcohol. All 55 animals were returned to the pasture to continue living after experimentation. All genomic DNA samples were extracted by using a QIAGEN DNeasy Blood & Tissue kit in accordance with the manufacturer’s protocol. Sequencing libraries were constructed with DNA extracts and a NEBNext ® Ltra DNA library preparation kit (Illumina ®, US). Sequencing was performed on an Illumina HiSeq × Ten platform (pair-end 150 bp). The sequencing data generated in this study were deposited in the NCBI SRA database (SRR10051499-SRR10551533 and SRR10613872-SRR10613891). In addition, we downloaded 166 non-intersex goat genome sequences from the NCBI SRA database. The detailed information of the 221 animals used in this study is shown in Supplementary Table 1 .

Read Filtering, Read Alignment, and Variant Calling

Raw sequencing reads were trimmed and filtered by using Trimmomatic (version 0. 36). We then mapped the clean pair-end reads to a goat ( Capra hircus ) reference genome (ARS1) by using BWA-MEM (version 0. 7. 13) with default parameters except that “-M” was enabled to flag shortened split hits as secondary data. We used Picard (version 2. 1. 1, http://broadinstitute. github. io/picard ) to remove potential PCR duplicates. Finally, the reads were locally realigned around indels with the IndelRealigner procedure in GATK (version 3. 7). We applied GATK to call variants and used the HaplotypeCaller algorithm in Genomic Variant Call Format (GVCF) mode. Variants were called individually for each animal, and one GVCF file that listed genotype likelihoods was generated per animal. Then, the variants were called from the GVCF files through joint genotyping analysis. We removed SNPs that were within the three base pairs of an indel by utilizing bcftools (version 1. 8). Biallelic SNPs were retained by applying a hard filter of QD < 2. 0, MQ < 40. 0, FS > 60. 0, SOR > 3. 0, MQRankSum < −12. 5, or ReadPosRankSum < −8. 0. We also used vcftools (version 0. 1. 14) to remove SNPs with a missing rating of more than 0. 1. The copy number variations (CNVs) with a silhouette score of <0. 65 and a MAF of <0. 05 were identified by using CNVcaller software ( Wang et al., 2017a ).

Genome-Wide Selective Sweep Analysis and Gene Annotation

Here, we carried out whole-genome selection signal analysis with two groups: (1) 35 intersex goats (case group, including 26 intersex Tangshan dairy goats [X/X] and 9 intersex Chinese goats [X/X] from southern China) vs. 186 non-intersex individuals (control group) and (2) 26 intersex goats (case group) vs. 20 non-intersex individuals (X/X, control group) from the Tangshan dairy goat population. For the SNP dataset, we calculated the pairwise fixation index ( F ST ) and π ratio (π intersex non−intersex ) with 40 kb sliding windows and 10 kb step size. Only windows passing the above two thresholds were retained. Candidate genes were subjected to functional enrichment with an online tool (KEGG, http://www. genome. jp/kegg/pathway. html ).

Additionally, we calculated V ST and F ST on the basis of absolute copy number (CN) to identify divergent CNV profiles between XX intersex and normal female goats. V ST was calculated by using following formula:

V S T = V t o t a l ( V p o p 1 × N p o p 1 + V p o p 2 × N p o p 2 )/ N t o t a l V t o t a l

, where V total is the total variance, V pop is the CN variance for each population, N pop is the sample size for each population, and N total is the total sample size.

Lastly, we calculated linkage disequilibrium by using Arlequin software version 3. 5. 1. 3 ( Excoffier and Lischer, 2010 ) with a permutation test (EM algorithm, permutation number = 100, 000).

PCR Amplification to Verify Structural Variant Genotypes

The primers ( Supplementary Table 2 ) of breakpoints based on the ~0. 48 Mb fragment (CHI1: 150, 334, 286–150, 818, 099) that was reversely inserted into the PIS deletion region (~10. 16 kb, CHI1: 129, 424, 780–129, 434, 939) were designed with the online tool Primer 3. 04 software ( http://bioinfo. ut. ee/primer3-0. 4. 0/ ) to identify the structural characteristics of the identified duplication variant in CHI1 in intersex goats. 2 × TransTaq ® High Fidelity (HiFi) PCR SuperMix II (TransGen Biotech, China) was used in PCR amplification. The qPCR reaction conditions consisted of an initial denaturation at 94°C for 5 min, followed by 35 cycles of denaturation at 94°C for 30 s, annealing at the locus-specific temperatures presented in Supplementary Table 2 for 30 s, and extension at 72°C for 120 s. Finally, an elongation step (final extension) was performed at 72°C for 7 min. Two-way Sanger sequencing was performed on an ABI 3730 sequencer platform (Life Technologies, US).

Three-Dimensional Genome Structure Comparison Between Intersex and Non-intersex Goats

Hi-C was performed on two individuals (Dazu black goat, China). One was intersex (hornless, PIS +/+), and the other was non-intersex (horned, female, PIS –/–). Sample processing (treatment of leukocytes from venous blood with the cell cross linker paraformaldehyde) and library construction were performed by using standard methods (four cutter restriction enzyme [MboI], Belton et al., 2012 ). Briefly as: (1) treat cells with paraformaldehyde (37% formaldehyde) to fix the conformation of DNA; (2) treat cross linked DNA with restriction enzymes (four cutter restriction enzyme, MboI) to produce sticky ends; (3) repair DNA ends with biotin labeling; (4) connect the DNA fragments by DNA ligase; (5) release the cross linked DNA state with 2. 5M Glycine; (6) purify the DNA by AMPure XP system (Beckman Coulter, Beverly, USA) and randomly break into 300~500 bp fragments; (7) construction of small DNA fragment library using NEB Next Ultra DNA Library Prep Kit (NEB, USA). After library construction, Qubit2. 0 was used for preliminary quantification. Then, the library was diluted to 1 ng /μL. Agilent 2100 was used to determine whether the insert size of the library met expectations. Q-PCR was used to quantify accurately the effective concentration of the library (> 2 nM), and sequencing was finally performed with an Illumina HiSeq × Ten PE150 TM platform. Sequencing data quality control, reference genome alignment (ARS1), interaction map construction, and loop structure analysis were performed with Juicer software ( Durand et al., 2016 ) with the standard parameters (Mbol restriction enzyme chunk size set at: 80000000 bp). Image visualization was performed by using the matplotlib package in the Python environment.

Results

A total of 16 462, 769 SNPs and 1, 058 CNVs were obtained from 221 samples. For the genome-wide selection of SNPs in all individuals (35 intersex vs. 186 non-intersex goats), we screened 258, 064 windows and estimated their F ST and π ratios ( Figure 1A , Supplementary Table 3 ). In total, we identified 40 windows in accordance with the intersection of the top 1% selective regions of both parameters (F ST and π ratio). These regions included 74 coding genes, which encompassed or were located up- and downstream within the 300 kb range of the window. Six of these genes were annotated to seven known signaling pathways ( Supplementary Table 4 ), including neuroactive ligand-receptor interaction ( P2RY13 and P2RY14 ), hippo signaling pathway-multiple species ( STK3 ), thyroid hormone signaling pathway ( MED12L ), hippo signaling pathway ( STK3 ), MTOR signaling pathway ( RRAGB ), protein processing in endoplasmic reticulum ( UBQLN2 ), and MAPK signaling pathway ( STK3 ).

FIGURE 1 Identification of a Goat Intersexuality-Associated Novel Variant Through Genome-Wide Resequencing and Hi-C Picture 1

Genome-wide selective sweep of goat intersexuality by using SNPs and CNVs.(A)Manhattan plot showing the SNP-based selection signals of 35 intersex goats compared with those of 186 non-intersex goats from a large geographically distributed population.(B)Manhattan plot showing the SNP-based selection signal of intersex goats within the Tangshan dairy goat population.(C)Manhattan plot showing the CNV-based selection signals of 35 intersex goats compared with those of 186 non-intersex goats from a large geographically distributed population.(D)Manhattan plot showing the CNV-based selection signals of intersex goats within the Tangshan dairy goat population.

Furthermore, we selected 46 Tangshan dairy goats from large samples, set up a scientific case–control analysis test, and identified the selected signal regions of intersexuality within Tangshan dairy goat populations to prevent the genetic background divergence of large sample populations from interfering with selection signal analysis ( Figure 1B , Supplementary Table 5 ). The results revealed that 50 windows were generated by the intersection of the top 1% selective regions of the F ST and π ratios and 79 genes, which covered or were located up- and downstream within the 300 kb range of the window. Only four genes were enriched in known pathways ( Supplementary Table 6 ). These genes included RBP2 (vitamin digestion and absorption), NCK1 ( ErbB signaling pathway, T cell receptor signaling pathway, and Axon guidance), IL20RB (Jak- STAT signaling pathway and cytokine–cytokine receptor interaction), and MRAS (tight junction, phospholipase D signaling pathway, proteoglycans in cancer, Rap1 signaling pathway, regulation of actin cytoskeleton, Ras and MAPK signaling pathways, and HTLV-I infection). Numerous consecutive windows in CHI1 (~129 to ~132 Mb) of intersex Tangshan dairy goats showed strong F ST and π ratio signals, and these windows covered the 11. 7 kb fragment deletion (PIS) that is widely recognized genomic signature of XX intersex goats ( Pailhoux et al., 2001 ).

The selective CNV-based sweep analysis of intersexuality with a large population and various genetic backgrounds of non-intersex goats (35 vs. 186) revealed that two CNVs had the highest signals with F ST (V1: F ST = 0. 834565, CHI1: 129, 424, 780–129, 434, 939; V2: F ST = 0. 830614, CHI1: 150, 334, 286–150, 818, 099) and V ST (V1–V ST = 0. 73290; V2– V ST = 0. 74050) ( Figure 1C , Supplementary Table 5 ). These two variants also carried the most prominent signal in the Tangshan dairy goat population (20 vs. 26; V1: F ST = 0. 60641, V ST = 0. 81321; V2: F ST = 0. 60641, V ST = 0. 83742; Figure 1D , Supplementary Table 7 ). The V1 variant was contained within the known intersex-related variant region (PIS deletion). Two genes, namely, MRPS22 (~120 kb distance) and COBB2 (~140 kb distance), were found upstream of V1, whereas the FOXL2 gene was found 340 kb downstream of V1. Other seven noncoding RNA ( LOC102190268, LOC108636915, LOC102185085, LOC108636375, LOC102190822, LOC102191084 , and LOC100861210 ) were discovered between FOXL2 and V1. Furthermore, KCNJ15 and ERG were encompassed by the V2 variant region, and ETS2 was located 230 kb downstream. No coding gene was found within the 500 kb upstream region of V2.

The length of the PIS deletion (V1) on CHI1 was ~10. 16 kb and was located from 129, 424, 780 bp to 129, 434, 939 bp ( Figures 2A, B ) as observed by using the IGV browser ( Thorvaldsdóttir et al., 2013 ). We found that the length of the V2 variant was ~0. 48 Mb and that this variant was distributed on CHI1 at 150, 334, 286–150, 818, 099 ( Figure 2C ). The different genotypes of the V1 and V2 variant regions could be clearly identified by comparing the read average coverage of each variant’s region with that of the whole genome ( Figures 2D, E ). On the basis of the genotypes of V1 and V2, we found that the homozygous deletion of V1 and the homozygous duplication of a ~0. 48 Mb region of V2 were always simultaneously present in all intersex goats ( Figure 2F ). Linkage disequilibrium analysis revealed significant linkage ( P < 0. 0001) between the V1 and V2 mutations in 221 goats.

FIGURE 2 Identification of a Goat Intersexuality-Associated Novel Variant Through Genome-Wide Resequencing and Hi-C Picture 2

Model detection and verification of novel transposition in intersex goats.(A)Two CNVs (V1 and V2) on Chromosome1 observed by IGV software.(B)Alignment and coverage of wide-genome short-reads from intersex goats as observed with IGV software. The CNV from 129, 424, 780 to 129, 434, 939 bp on Chromosome 1 (V1, B) manifested as a deletion.(C)Alignment and coverage of wide-genome short-reads from intersex goats as observed by IGV software. The CNV from 150, 334, 286 to 150, 818, 099 bp on Chromosome 1 manifested as a duplication (V2).(D)Genomic coverage of different genotypes of V1 variant reads (Chromosome 1, 129. 40–129. 45 Mb).(E)Genomic coverage of different genotype of V2 variant reads (Chromosome 1, 150. 20–150. 90 Mb).(F)Two CNV variants (V1 and V2) associated with intersexuality had the same frequency in the population.(G)Schematic of the Chromosome 1 PIS transposition model and location map of primers for PCR verification.(H)Gel electrophoresis verification of PCR results.(I, J)Sanger sequencing results of sequences amplified with primers 4 and 5.

In the heterozygous and extra duplication homozygous individuals with the V2 mutation, a considerable number of reads were split-mapped simultaneously to the outer boundary of V1 and the inner boundary of V2 ( Supplementary Figure 1 ). We verified the true boundary breakpoints of the two variant regions through PCR amplification and Sanger sequencing ( Figures 2G–J ). Therefore, the precise PIS genome structure was doubly confirmed as an inverted duplication of the ~0. 48 Mb segment that had inserted into the 10. 16 kb PIS deletion.

An average of > 250 Gb (~85×) of genome coverage data were obtained from two individuals and used to construct a 3D genome high-resolution interaction map. Firstly, we gathered total 851, 483, 595 and 808, 343, 466 reads in case and control individual, respectively. Secondly, according to the mapping results of case and control dataset, there are 300, 895, 514 (35. 34%) and 359, 770, 248 (44. 51%) normal paired reads, 402, 929, 844 (47. 32%) and 347, 404, 121 (42. 98%) chimeric paired reads, 344, 755, 346 and 117, 458, 553 PCR duplicates reads, 156, 341, 608 and 244, 912, 442 intra-chromosomal interaction reads, 71, 982, 908 and 173, 316, 912 short-distance interaction sequence with interaction distance less than 20 kb (<20 kb), 84, 358, 664 and 149, 510, 017 long-range interaction sequences with distance larger than 20 kb (> 20 kb), respectively.

The heat map of both sets of data with 80 kb resolution revealed a potential intrachromosomal rearrangement site ( Figures 3A, B ), which was initially identified in CHI1 of an intersex individual (PIS–/–). This finding was consistent with the physical location of V2, which was absent from non-intersex goats (PIS–/–). We used a resolution of 10 kb to identify effective breakpoints ( Figures 3C, D ). A small but sharp contact peak suggested that a new intrachromosomal rearrangement occurred in CHI1 of the intersex goats ( Figures 3E, F ). We identified four private consequent loop regions in CHI1 of the intersex goats ( Supplementary Table 8 ) and compared these regions with those in non-intersex individuals ( Supplementary Table 9 ). These loop regions were densely clustered in the ~20 kb downstream regions of the FOXL2 gene, which overlapped with LOC102191651 and LOC108636917 ( Figure 3G ).

FIGURE 3 Identification of a Goat Intersexuality-Associated Novel Variant Through Genome-Wide Resequencing and Hi-C Picture 3

Hi-C analysis results.(A)Window interaction matrix on CHI1 from 128. 43 to 151. 08 Mb in intersexuality individuals with 80 Kb resolution.(B)Window interaction matrix on CHI1 from 128. 43 to 151. 08 Mb in normal control individuals with 80 Kb resolution.(C–F)Interaction matrix and comparison of the 129. 1–129. 7 and 150. 3–150. 9 Mb regions on CHI1 between intersexuality and normal individuals with 10 Kb resolution.(G)Analysis of loop conformation in the 129. 4–130. 2 Mb region of CHI1 between intersexuality and normal individuals at 5 kb resolution revealed that the intersexuality case had a particular loop conformation in the 129. 80–129. 88 Mb region at which the sexual development-associated gene FOXL2 was located 30 Kb upstream.

Discussion

SNP-based genome-wide selection signal analysis revealed numerous sharp signals in 35 intersex goats and 186 control samples. Within the top 1% selection window, a series of genes were identified and found to be deeply involved in animal reproduction and multiple developmental processes. For example, STK3 is a key molecule that connects the downstream signaling pathway of estrogen and the Hippo signaling pathway; it also regulates the dynamic development of the uterine epithelium during the estrous cycle through the signal transduction of uterine epithelial cells ( Moon et al., 2019 ). STK3 was annotated to other signaling pathways, such as the MAPK ( Bogani et al., 2009 ; Warr et al., 2016 ) and Hippo signaling pathways, that also play an important role in gonadal development and sex determination ( Frum et al., 2018 ; Devos et al., 2020 ). Although MED12L has been verified to be associated with fetal mental retardation in human ( Nizon et al., 2019 ), it is also involved in reproductive development ( Sayem et al., 2017 ; Das and Kumar, 2018 ; Ulloa-Aguirre et al., 2018 ). In addition, the RRAGB gene is enriched in the mTOR signaling pathway, which is widely involved in gonadal development ( Bajwa et al., 2017 ; Correia et al., 2020 ). Therefore, our findings suggested that numerous molecular mechanisms underlying development and the physiological maintenance of intersexual characteristics await further excavation.

Notably, the different genetic backgrounds of large samples can cause many false-positive genes, and chromosomal regions may thus be incorrectly identified. Given the inconsistent ratios between the numbers of Tangshan dairy goats in the intersex and control groups, a gene with considerable breed specificity caused interference. We performed a strict case–control experiment on the Tangshan dairy goat population to prevent the interference of specific population backgrounds and identified a series of interesting genes. The highest continuous selection signal was observed in CHI1. These signals covered the areas of a previously reported PIS deletion ( Pailhoux et al., 2001 ) and six coding genes ( MRAS, NMNAT3, ARMC8, DBR1, LOC108636376 , and SCLC35G2 ); SOX14 and MRPS22 in the upstream region; and LOCl02190268 and IL20RB in the downstream region.

Some genes with cellular biological importance were identified. For example, NMNAT3 maintains cell differentiation by maintaining mitochondrial content ( Son et al., 2016 ; Yu et al., 2020 ). ARMC8 is involved in the adherence of regulatory cells to cells and is associated with cell differentiation in ovarian cancer tumors ( Jiang et al., 2015 ; Gul et al., 2019 ). SOX14 is associated with apoptosis in cancer cells in the sexual reproductive system ( Stanisavljevic et al., 2017 ) and is a crucial determinant of allergy development in Drosophila ( Ritter and Beckstead, 2010 ). Interestingly, the conserved region of the MRPS22 gene is a long-range enhancer and regulates the expression of FOXL2 through an unclear advanced cis -regulatory effect of chromatin structure in humans and rats ( Crisponi et al., 2004 ). Thus, our case–control experiment involving a population with a single genetic background enabled us to screen out many false-positive signals and identify a series of credible candidate genes. The results of this experiment provided insight into the molecular genetic mechanism of intersexuality-related physiology.

The recent identification of the gene transcription profiles of intersex and normal goat gonads through the use of RNA-Seq technology suggests that a large number of differentially expressed genes may be involved in the regulation of sex determination and differentiation in intersex goats ( Yang et al., 2020 ). This result reminds us that many potential molecular mechanisms under the goat sexual reversal phenotype remain unclear.

Our CNV-based analysis results showed that equally strong signals were generated in V1 and V2 in the large sample with different genetic backgrounds and the Tangshan dairy goat population in the case–control analysis. These signals were recognized as the 10. 16 kb PIS deletion and the ~0. 4838 Mb duplicated segment located ~20. 9 Mb further downstream of the PIS deletion and ~150 Mb on CHI1. In addition, as expected, the highly complex structure was identified as the additional 0. 4838 Mb-sized duplicated segment that was inversely inserted at the breakpoint of the 10. 16 kb deletion. Our findings were consistent with the recent research results from a team in Germany that used long-read whole-genome sequencing ( Simon et al., 2020 ). Although we utilized short-read sequencing technology, the large sample size and classic case–control experimental design still achieved the same effect. Our study confirmed that the XX intersex goats from the hornless goat population in China share the same PIS genome variant structure with European goats.

In accordance with previous studies that identified the segment size and polymorphism in PIS deletion ( Li et al., 2011 ). We believe that some loss in the ASR1 assembly occurred on the last 180 bp section, while was not lost on the previous 11. 7 kb PIS deletion sequence (GenBank No. AF404302) investigated by Simon et al. (2020). It was adjacent to the PISRT1 gene with the closest physical distance. However, previous studies have shown that PISRT1 does not participate in the expression of FOXL2 and the determination/differentiation of the gonads. For example, the overexpression of PISRT1 in PIS–/– fetuses does not affect FOXL2 expression levels and gonadal development ( Boulanger et al., 2008 ).

The duplicated segment contained the KCNJ15 and ERG genes. The extra copies of these two genes have an essential role in horn and gonadal development. KCNJ15 is known to participate in insulin secretion ( Okamoto et al., 2012 ), nervous system diseases ( Zhou et al., 2018 ), gastric acid secretion ( Yuan et al., 2015 ), kidney cancer ( Liu et al., 2019 ), and esophageal squamous cell carcinoma ( Nakamura et al., 2020 ). It is also involved in gastric acid secretion and regulation ( He et al., 2011 ), and the relationship between gastric acid secretion and the effects of sex hormones was verified decades ago ( Adeniyi, 1991 ). The high expression of KCNJ15 in follicle-associated epithelium suggests that KCNJ15 may be involved in the functional development of the ovary ( Kobayashi et al., 2012 ) and implicates this gene in female gonadal development. Furthermore, a large number of studies have shown that ERG is not only an oncogene that is related to a variety of cancers ( Wang et al., 2017b ; Zhang et al., 2020 ), it also participates in the embryonic developmental processes, including bone development ( Iwamoto et al., 2000 ), of a variety of organisms ( Furlan et al., 2005 ; Nikolova-Krstevski et al., 2009 ). This participation indicates that the ERG gene may be related to horn and embryonic development.

Furthermore, Hi-C technology was used to study DNA replication, transcription regulation, and DNA damage repair and contact between chromosomal loci ( Cremer and Cremer, 2001 ; de Wit and De Laat, 2012 ; Maass et al., 2018 ). Currently, this topic is heavily explored in genomic research, and numerous studies on technical method optimization have been performed ( Lin et al., 2018 ; Yardimci et al., 2019 ; Janaratne et al., 2020 ).

Intrachromosomal rearrangement or palindrome duplication is associated with various processes of phenotypic determination and development ( Carbonell-Bejerano et al., 2017 ; Yin et al., 2017 ; Mendoza et al., 2020 ). We performed the loop analysis of the 3D genomes to further investigate the special chromosomal spatial structures resulting from the identified intrachromosomal rearrangement. We found several unique loop structures in CHI1 of homozygous PIS intersex goats but not in that of non-intersex individuals. Many intrachromosomal rearrangement structures can alter gene expression levels within and in areas adjacent to a gene region by altering chromosomal structure ( Demura et al., 2007 ; Suzuki et al., 2020 ).

Substantial evidence indicates that many of the observed loops are related to gene regulation and serve as anchors and promoters ( Ahmadiyeh et al., 2010 ; Hoffman et al., 2013 ; Rao et al., 2014 ). The loops that we identified in this study were consistent and clustered near the FOXL2 gene. Numerous pieces of evidence have shown that silencing FOXL2 expression directly affects ovarian development and oogenesis in fish ( Fan et al., 2019 ), mice, and humans ( Uda et al., 2004 ; Thanatsis et al., 2019 ). Specifically, the elimination of FOXL2 expression is sufficient to induce female-to-male reversal in XX goats ( Pannetier et al., 2012 ; Boulanger et al., 2014 ). Therefore, although the regulatory relationship between this newly discovered intrachromosomal rearrangement and FOXL2 expression cannot be evaluated thus far, the change in spatial chromosome 3D structure in the adjacent region of FOXL2 was evident. Whether these loop structures affect FOXL2 expression and cause intersexuality by inhibiting cis -acting elements or switching trans -acting elements should be evaluated through in-depth molecular biology research.

In addition, we found that two genes were located within the loop region: one was trafficking protein particle complex subunit 1 pseudogene (LOC102191651), and the other was an uncharacterized non-coding RNA (LOC108636917). Additional evidence regarding the further functions of these both genes remains lacking. Therefore, we cannot conclude that these loops/two genes participate in gonadal development. However, an interesting gene, PIK3CB , that was located further downstream of LOC102191651 and LOC108636917 attracted our attention. Numerous studies have shown that PIK3CB plays an important role in the development and physiological function of the ovary ( Zheng et al., 2012 ; Li et al., 2013 ; Nteeba et al., 2017 ). However, supernumerary data suggesting that this gene is responsible for the occurrence and maintenance of the intersexual phenotype are unavailable. Therefore, whether the novel loop region containing both genes affects PIK3CB expression and whether PIK3CB is a new essential factor that is sufficient for causing female-to-male sex reversal in XX goats need to be evaluated.

Conclusions

We performed the genome-wide selective sweep of intersex goats with wide-genome next-generation sequencing. We doubly verified that the structural variant of caprine PIS structure, a 0. 48 Mb duplicated fragment located ~20 Mb downstream of the PIS region that was reversely inserted into the PIS deletion, was sufficient as a broad-spectrum clinical diagnostic marker of XX intersex goats from Europe and China. The existence of several private dense loop structures in the region adjacent to FOXL2 of intersex XX goats but not in that of non-intersex individuals suggested that intrachromosomal rearrangement might affect the expression of FOXL2 or other neighboring novel candidate genes. This effect needs to be further evaluated. This study supported a precise genomic feature of PIS phenotype in intersex goats from Europe and China and provided new insights for future research on the molecular genetic mechanism underlying female-to-male sex reversal in goats.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found at: https://www. ncbi. nlm. nih. gov/genbank/ , SRR10051499-SRR10551533; https://www. ncbi. nlm. nih. gov/genbank/ , SRR10613872-SRR10613891.

Ethics Statement

The animal study was reviewed and approved by Animal Experiments of the Southwest University (No. [2007] 3) Southwest University, Chongqing 400716, China.

Author Contributions

G-XE, YJ, Y-FH, Y-JZ, Y-HM, J-HZ, Q-HH, M-XC, and X-LL conceived and designed the experiments. D-KZ, Z-QZ, B-GY, H-JG, Y-MH, and X-HD analyzed the data. G-XE analyzed the data and wrote the paper. YJ, X-PJ, R-SN, Y-GH, YZ, R-YZ, W-HN, and L-HL supported the samples and displayed lab work. G-XE, Y-FH, and YJ provided funding. All authors read and approved the manuscript.

Funding

This work was supported by the National Key R&D Program of China (No. 2018YFD0502000), Characteristic Germplasm Resources Population Selection and Innovation on Mutton Sheep and Goats (No. 2015BAD03B05), National Natural Science Foundation of China (No. 31172195), Chongqing Research Program of Basic Research and Frontier Technology (cstc2018jcyjAX0153), and Fundamental Research Funds for the Central Universities (XDJK2018B014). Chongqing Science and Technology Innovation Special Project (No. cstc2019jscx-gksbX0135). We are grateful to thank Prof. Dr. Gang Cao from Huazhong Agricultural University (Wuhan, China) for assisting in part of bioinformatics analysis of Hi-C datasets.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www. frontiersin. org/articles/10. 3389/fgene. 2020. 616743/full#supplementary-material

Abbreviations

PIS, Polled intersex syndrome; CNV, Copy number variant; FST, Pairwise fixation index; SNP, Single nucleotide polymorphism; Hi-C, High-through chromosome conformation capture; CHI, Capra hircus chromosome; FOXL2 , Forkhead box L2; KCNJ15 , Potassium inwardly rectifying channel subfamily J member 15; PIK3CB , Phosphatidylinositol-4, 5-bisphosphate 3-kinase catalytic subunit beta.

References

Adeniyi, K. O. (1991). Gastric acid secretion and parietal cell mass: effect of sex hormones. Gastroenterology 101, 66–69. doi: 10. 1016/0016-5085(91)90460-3

Ahmadiyeh, N., Pomerantz, M. M., Grisanzio, C., Herman, P., Jia, L., Almendro, V., et al. (2010). 8q24 prostate, breast, and colon cancer risk loci show tissue-specific long-range interaction with MYC. Proc. Natl. Acad. Sci. U. S. A . 107, 9742–9746. doi: 10. 1073/pnas. 0910668107

Bagheri-Fam, S., Bird, A. D., Zhao, L., Ryan, J. M., Yong, M., Wilhelm, D., et al. (2017). Testis determination requires a specific FGFR2 isoform to repress FOXL2. Endocrinology 158, 3832–3843. doi: 10. 1210/en. 2017-00674

Bajwa, P., Nielsen, S., Lombard, J. M., Rassam, L., Nahar, P., Rueda, B. R., et al. (2017). Overactive mTOR signaling leads to endometrial hyperplasia in aged women and mice. Oncotarget 8: 7265. doi: 10. 18632/oncotarget. 13919

Baron, D., Batista, F., Chaffaux, S., Cocquet, J., Cotinot, C., Cribiu, E., et al. (2005). Foxl2 gene and the development of the ovary: a story about goat, mouse, fish and woman. Reprod. Nutr. Dev. 45, 377–382. doi: 10. 1051/rnd: 2005028

Belton, J. M., McCord, R. P., Gibcus, J. H., Naumova, N., Zhan, Y., and Dekker, J. (2012). Hi-C: a comprehensive technique to capture the conformation of genomes. Methods 58, 268–276. doi: 10. 1016/j. ymeth. 2012. 05. 001

Bogani, D., Siggers, P., Brixey, R., Warr, N., Beddow, S., Edwards, J., et al. (2009). Loss of mitogen-activated protein kinase kinase kinase 4 (MAP3K4) reveals a requirement for MAPK signalling in mouse sex determination. PLoS Biol 7: e1000196. doi: 10. 1371/journal. pbio. 1000196

Bosu, W., and Basrur, P. K. (1984). Morphological and hormonal features of an ovine and a caprine intersex. Can. J. Comp. Med. 48: 402.

Boulanger, L., Kocer, A., Daniel, N., Pannetier, M., Chesné, P., Heyman, Y., et al. (2008). Attempt to rescue sex-reversal by transgenic expression of the PISRT1 gene in XX PIS–/–goats. Sexual Development 2, 142–151. doi: 10. 1159/000143432

Boulanger, L., Pannetier, M., Gall, L., Allais-Bonnet, A., Elzaiat, M., Le Bourhis, D., et al. (2014). FOXL2 is a female sex-determining gene in the goat. Curr. Biol. 24, 404–408. doi: 10. 1016/j. cub. 2013. 12. 039

Carbonell-Bejerano, P., Royo, C., Torres-Pérez, R., Grimplet, J., Fernandez, L., Franco-Zorrilla, J. M., et al. (2017). Catastrophic unbalanced genome rearrangements cause somatic loss of berry color in grapevine. Plant Physiol. 175, 786–801. doi: 10. 1104/pp. 17. 00715

Chen, Y., Chen, J., and Jiang, H. Z. (2010). Research on the characteristics of intersex goats and it’s related research progress. J. China Anim. Husbandry Vet. 37, 125–129 (in Chinese).

Correia, B., Sousa, M. I., and Ramalho-Santos, J. (2020). The mTOR pathway in reproduction: from gonadal function to developmental coordination. Reproduction 159, R173–R188. doi: 10. 1530/REP-19-0057

Cremer, T., and Cremer, C. (2001). Chromosome territories, nuclear architecture and gene regulation in mammalian cells. Nat. Rev. Genet. 2, 292–301. doi: 10. 1038/35066075

Crisponi, L., Uda, M., Deiana, M., Loi, A., Nagaraja, R., Chiappe, F., et al. (2004). FOXL2 inactivation by a translocation 171 kb away: analysis of 500 kb of chromosome 3 for candidate long-range regulatory sequences. Genomics 83, 757–764. doi: 10. 1016/j. ygeno. 2003. 11. 010

Das, N., and Kumar, T. R. (2018). Molecular regulation of follicle-stimulating hormone synthesis, secretion and action. J. Mol. Endocrinol. 60, R131–R155. doi: 10. 1530/JME-17-0308

de Wit, E., and De Laat, W. (2012). A decade of 3C technologies: insights into nuclear organization. Genes Dev. 26, 11–24. doi: 10. 1101/gad. 179804. 111

Demura, M., Martin, R. M., Shozu, M., Sebastian, S., Takayama, K., Hsu, W.-T., et al. (2007). Regional rearrangements in chromosome 15q21 cause formation of cryptic promoters for the CYP19 (aromatase) gene. Hum. Mol. Genet. 16, 2529–2541. doi: 10. 1093/hmg/ddm145

Devos, M., Grosbois, J., and Demeestere, I. (2020). Interaction between PI3K/AKT and Hippo pathways during in vitro follicular activation and response to fragmentation and chemotherapy exposure using a mouse immature ovary model. Biol. Reprod. 102, 717–729. doi: 10. 1093/biolre/ioz215

Durand, N. C., Shamim, M. S., Machol, I., Rao, S. S., Huntley, M. H., Lander, E. S., et al. (2016). Juicer provides a one-click system for analyzing loop-resolution Hi-C experiments. Cell Syst. 3, 95–98. doi: 10. 1016/j. cels. 2016. 07. 002

Eaton, O. N. (1945). The relation between polled and hermaphroditic characters in dairy goat. Genetics 30, 51–61.

Elzaiat, M., Todeschini, A. L., Caburet, S., and Veitia, R. (2017). The genetic make-up of ovarian development and function: the focus on the transcription factor FOXL2. Clin. Genet. 91, 173–182. doi: 10. 1111/cge. 12862

Excoffier, L., and Lischer, H. E. (2010). Arlequin suite ver 3. 5: a new series of programs to perform population genetics analyses under Linux and Windows. Mol. Ecol. Resour. 10, 564–567. doi: 10. 1111/j. 1755-0998. 2010. 02847. x

Fan, Z., Zou, Y., Liang, D., Tan, X., Jiao, S., Wu, Z., et al. (2019). Roles of forkhead box protein L2 (foxl2) during gonad differentiation and maintenance in a fish, the olive flounder (Paralichthys olivaceus). Reprod. Fertil. Dev. 31, 1742–1752. doi: 10. 1071/RD18233

Frum, T., Murphy, T. M., and Ralston, A. (2018). HIPPO signaling resolves embryonic cell fate conflicts during establishment of pluripotency in vivo . Elife 7: e42298. doi: 10. 7554/eLife. 42298. 020

Furlan, F., Guasti, L., Avossa, D., Becchetti, A., Cilia, E., Ballerini, L., et al. (2005). Interneurons transiently express the ERG K+ channels during development of mouse spinal networks in vitro . Neuroscience 135, 1179–1192. doi: 10. 1016/j. neuroscience. 2005. 06. 040

Gul, I. S., Hulpiau, P., Sanders, E., Van Roy, F., and van Hengel, J. (2019). Armc8 is an evolutionarily conserved armadillo protein involved in cell–cell adhesion complexes through multiple molecular interactions. Biosci. Rep. 39: BSR20180604. doi: 10. 1042/BSR20180604

He, W., Liu, W., Chew, C. S., Baker, S. S., Baker, R. D., Forte, J. G., et al. (2011). Acid secretion-associated translocation of KCNJ15 in gastric parietal cells. Am. J. Physiol. Gastrointest. Liver Physiol. 301, G591–G600. doi: 10. 1152/ajpgi. 00460. 2010

Hoffman, M. M., Ernst, J., Wilder, S. P., Kundaje, A., Harris, R. S., Libbrecht, M., et al. (2013). Integrative annotation of chromatin elements from ENCODE data. Nucleic Acids Res. 41, 827–841. doi: 10. 1093/nar/gks1284

Iwamoto, M., Higuchi, Y., Koyama, E., Enomoto-Iwamoto, M., Kurisu, K., Yeh, H., et al. (2000). Transcription factor ERG variants and functional diversification of chondrocytes during limb long bone development. J. Cell Biol. 150, 27–40. doi: 10. 1083/jcb. 150. 1. 27

Janaratne, T., Wang, X. C., Becker, C., Zhao, Y., Leanna, R., and Pritts, W. (2020). Hydrophobic interaction chromatography (HIC) method development and characterization of resolved drug-load variants in site-specifically conjugated pyrrolobenzodiazepine dimer-based antibody drug conjugates (PBD-ADCs). J. Pharm. Biomed. Anal. 179: 113027. doi: 10. 1016/j. jpba. 2019. 113027

Jiang, G., Yang, D., Wang, L., Zhang, X., Xu, H., Miao, Y., et al. (2015). A novel biomarker ARMc8 promotes the malignant progression of ovarian cancer. Hum. Pathol. 46, 1471–1479. doi: 10. 1016/j. humpath. 2015. 06. 004

Kijas, J. W., Ortiz, J. S., McCulloch, R., James, A., Brice, B., Swain, B., et al. (2013). Genetic diversity and investigation of polledness in divergent goat populations using 52 088 SNPs. Anim. Genet. 44, 325–335. doi: 10. 1111/age. 12011

Kobayashi, A., Donaldson, D. S., Kanaya, T., Fukuda, S., Baillie, J. K., Freeman, T. C., et al. (2012). Identification of novel genes selectively expressed in the follicle-associated epithelium from the meta-analysis of transcriptomics data from multiple mouse cell and tissue populations. DNA Res. 19, 407–422. doi: 10. 1093/dnares/dss022

Li, Q., He, H., Zhang, Y.-L., Li, X.-M., Guo, X., Huo, R., et al. (2013). Phosphoinositide 3-kinase p110δ mediates estrogen-and FSH-stimulated ovarian follicle growth. Mol. Endocrinol. 27, 1468–1482. doi: 10. 1210/me. 2013-1082

Li, X., Zhang, J., Zhou, R., Li, L., and Zheng, G. (2011). Special variations within 11. 7 kb fragment in goat polled intersex syndrome. Afr. J. Biotechnol. 10, 6695–6699. doi: 10. 5897/AJB11. 001

Lin, D., Hong, P., Zhang, S., Xu, W., Jamal, M., Yan, K., et al. (2018). Digestion-ligation-only Hi-C is an efficient and cost-effective method for chromosome conformation capture. Nat. Genet. 50, 754–763. doi: 10. 1038/s41588-018-0111-2

Liu, Y., Wang, H., Ni, B., Zhang, J., Li, S., Huang, Y., et al. (2019). Loss of KCNJ15 expression promotes malignant phenotypes and correlates with poor prognosis in renal carcinoma. Cancer Manag. Res. 11: 1211. doi: 10. 2147/CMAR. S184368

Maass, P. G., Barutcu, A. R., Weiner, C. L., and Rinn, J. L. (2018). Inter-chromosomal contact properties in live-cell imaging and in Hi-C. Mol. Cell 69, 1039–1045. e3. doi: 10. 1016/j. molcel. 2018. 02. 007

Mendoza, M. N., Schalnus, S. A., Thomson, B., Bellone, R. R., Juras, R., and Raudsepp, T. (2020). Novel complex unbalanced dicentric X-autosome rearrangement in a thoroughbred mare with a mild effect on the phenotype. Cytogen. Genome Res . 160, 597–609. doi: 10. 1159/000511236

Moon, S., Lee, O.-H., Lee, S., Lee, J., Park, H., Park, M., et al. (2019). STK3/4 expression is regulated in uterine endometrial cells during the estrous cycle. Cells 8: 1643. doi: 10. 3390/cells8121643

Nakamura, S., Kanda, M., Koike, M., Shimizu, D., Umeda, S., Hattori, N., et al. (2020). KCNJ15 expression and malignant behavior of esophageal squamous cell carcinoma. Ann. Surg. Oncol . 27, 2559–2568. doi: 10. 1245/s10434-019-08189-8

Nikolova-Krstevski, V., Yuan, L., Le Bras, A., Vijayaraj, P., Kondo, M., Gebauer, I., et al. (2009). ERG is required for the differentiation of embryonic stem cells along the endothelial lineage. BMC Dev. Biol 9: 72. doi: 10. 1186/1471-213X-9-72

Nizon, M., Laugel, V., Flanigan, K. M., Pastore, M., Waldrop, M. A., Rosenfeld, J. A., et al. (2019). Variants in MED12L, encoding a subunit of the mediator kinase module, are responsible for intellectual disability associated with transcriptional defect. Genet. Med . 21, 2713–2722. doi: 10. 1038/s41436-019-0557-3

Nteeba, J., Ganesan, S., Madden, J. A., Dickson, M. J., and Keating, A. F. (2017). Progressive obesity alters ovarian insulin, phosphatidylinositol-3 kinase, and chemical metabolism signaling pathways and potentiates ovotoxicity induced by phosphoramide mustard in mice. Biol. Reprod. 96, 478–490. doi: 10. 1095/biolreprod. 116. 143818

Okamoto, K., Iwasaki, N., Doi, K., Noiri, E., Iwamoto, Y., Uchigata, Y., et al. (2012). Inhibition of glucose-stimulated insulin secretion by KCNJ15, a newly identified susceptibility gene for type 2 diabetes. Diabetes 61, 1734–1741. doi: 10. 2337/db11-1201

Pailhoux, E., Vigier, B., Chaffaux, S., Servel, N., Taourit, S., Furet, J.-P., et al. (2001). A 11. 7-kb deletion triggers intersexuality and polledness in goats. Nat. Genet. 29, 453–458. doi: 10. 1038/ng769

Pannetier, M., Chassot, A.-A., Chaboissier, M.-C., and Pailhoux, E. (2016). Involvement of FOXL2 and RSPO1 in ovarian determination, development, and maintenance in mammals. Sex. Dev. 10, 167–184. doi: 10. 1159/000448667

Pannetier, M., Elzaiat, M., Thépot, D., and Pailhoux, E. (2012). Telling the story of XX sex reversal in the goat: highlighting the sex-crossroad in domestic mammals. Sex. Dev. 6, 33–45. doi: 10. 1159/000334056

Ramadan, R., and El Hassan, A. M. (1988). Intersexuality in goats. N. Z. Vet. J. 36, 120–124. doi: 10. 1080/00480169. 1988. 35505

Ramadan, R., Gameel, A., Dafalla, E., and Galil, A. (1991). Bilateral scrotal hysterocele in an intersex goat. J. Vet. Med. Ser. A 38, 441–444. doi: 10. 1111/j. 1439-0442. 1991. tb01033. x

Rao, S. S., Huntley, M. H., Durand, N. C., Stamenova, E. K., Bochkov, I. D., Robinson, J. T., et al. (2014). A 3D map of the human genome at kilobase resolution reveals principles of chromatin looping. Cell 159, 1665–1680. doi: 10. 1016/j. cell. 2014. 11. 021

Ritter, A. R., and Beckstead, R. B. (2010). Sox14 is required for transcriptional and developmental responses to 20-hydroxyecdysone at the onset of Drosophila metamorphosis. Dev. Dyn. 239, 2685–2694. doi: 10. 1002/dvdy. 22407

Sayem, A. S. M., Giribabu, N., Muniandy, S., and Salleh, N. (2017). Effects of thyroxine on expression of proteins related to thyroid hormone functions (TR-α, TR-β, RXR and ERK1/2) in uterus during peri-implantation period. Biomed. Pharmacother. 96, 1016–1021. doi: 10. 1016/j. biopha. 2017. 11. 128

Schibler, L., Cribiu, E. P., Oustry-Vaiman, A., Furet, J.-P., and Vaiman, D. (2000). Fine mapping suggests that the goat Polled Intersex Syndrome and the human Blepharophimosis Ptosis Epicanthus Syndrome map to a 100-kb homologous region. Genome Res. 10, 311–318. doi: 10. 1101/gr. 10. 3. 311

Simon, R., Lischer, H., Pieńkowska-Schelling, A., Keller, I., Häfliger, I. M., Letko, A., et al. (2020). New genomic features of the polled intersex syndrome variant in goats unraveled by long-read whole-genome sequencing. Anim. Genet. 51, 439–448. doi: 10. 1111/age. 12918

Son, M. J., Kwon, Y., Son, T., and Cho, Y. S. (2016). Restoration of mitochondrial NAD+ levels delays stem cell senescence and facilitates reprogramming of aged somatic cells. Stem Cells 34, 2840–2851. doi: 10. 1002/stem. 2460

Song, C. F. L. Y., Zhang, D. X., Jiao, L., Wang, G. M., and Liu, X. J. (2015). Anatomical observation of reproductive organs in intersex goats. Heilongjiang Anim. Sci. Vet. Med . 16: 81. (in Chinese).

Stanisavljevic, D., Petrovic, I., Vukovic, V., Schwirtlich, M., Gredic, M., Stevanovic, M., et al. (2017). SOX14 activates the p53 signaling pathway and induces apoptosis in a cervical carcinoma cell line. PLoS ONE 12: e0184686. doi: 10. 1371/journal. pone. 0184686

Suzuki, M., Katayama, S., and Yamamoto, M. (2020). Two effects of GATA2 enhancer repositioning by 3q chromosomal rearrangements. IUBMB Life 72, 159–169. doi: 10. 1002/iub. 2191

Tao, W., Chen, J., Tan, D., Yang, J., Sun, L., Wei, J., et al. (2018). Transcriptome display during tilapia sex determination and differentiation as revealed by RNA-Seq analysis. BMC Genomics 19: 363. doi: 10. 1186/s12864-018-4756-0

Thanatsis, N., Kaponis, A., Koika, V., Georgopoulos, N. A., and Decavalas, G. O. (2019). Reduced Foxo3a, FoxL2, and p27 mRNA expression in human ovarian tissue in premature ovarian insufficiency. Hormones 18, 409–415. doi: 10. 1007/s42000-019-00134-4

Thorvaldsdóttir, H., Robinson, J. T., and Mesirov, J. P. (2013). Integrative Genomics Viewer (IGV): high-performance genomics data visualization and exploration. Brief. Bioinformatics 14, 178–192. doi: 10. 1093/bib/bbs017

Uda, M., Ottolenghi, C., Crisponi, L., Garcia, J. E., Deiana, M., Kimber, W., et al. (2004). Foxl2 disruption causes mouse ovarian failure by pervasive blockage of follicle development. Hum. Mol. Genet. 13, 1171–1181. doi: 10. 1093/hmg/ddh124

Ulloa-Aguirre, A., Zariñán, T., Jardón-Valadez, E., Gutiérrez-Sagal, R., and Dias, J. A. (2018). Structure-function relationships of the follicle-stimulating hormone receptor. Front. Endocrinol. 9: 707. doi: 10. 3389/fendo. 2018. 00707

Vaiman, D., Koutita, O., Oustry, A., Elsen, J.-M., Manfredi, E., Fellous, M., et al. (1996). Genetic mapping of the autosomal region involved in XX sex-reversal and horn development in goats. Mamm. Genome 7, 133–137. doi: 10. 1007/s003359900033

Vaiman, D., Schibler, L., Oustry-Vaiman, A., Pailhoux, E., Goldammer, T., Stevanovic, M., et al. (1999). High-resolution human/goat comparative map of the goat polled/intersex syndrome (PIS): the human homologue is contained in a human YAC from HSA3q23. Genomics 56, 31–39. doi: 10. 1006/geno. 1998. 5691

Wang, X., Qiao, Y., Asangani, I. A., Ateeq, B., Poliakov, A., Cieślik, M., et al. (2017a). Development of peptidomimetic inhibitors of the ERG gene fusion product in prostate cancer. Cancer Cell 31, 532–548. e7. doi: 10. 1016/j. ccell. 2017. 02. 017

Wang, X., Zheng, Z., Cai, Y., Chen, T., Li, C., Fu, W., et al. (2017b). CNVcaller: highly efficient and widely applicable software for detecting copy number variations in large populations. GigaScience 6: gix115. doi: 10. 1093/gigascience/gix115

Wang, Y. Q., and Zhang, X. (1993). Study on cytogenetics of intersex livestock. J. Zhengzhou Pastoral Sci. 13, 20–23.

Warr, N., Siggers, P., Carré G.-A, Wells, S., and Greenfield, A. (2016). Genetic analyses reveal functions for MAP2K3 and MAP2K6 in mouse testis determination. Biol. Reprod. 94, 1–7. doi: 10. 1095/biolreprod. 115. 138057

Yang, B., Jia, L., Zhao, D., Meng, L., Liu, X., Zhang, Y., et al. (2012). Detection and application of PIS genetic deficiency gene in dairy goat. Hereditas 34, 895–900. doi: 10. 3724/SP. J. 1005. 2012. 00895

Yang, S., Han, H., Li, J., Zhang, Y., Zhao, J., Wei, H., et al. (2020). Transcriptomic analysis of gene expression in normal goat ovary and intersex goat gonad. Reprod. Domestic Anim . doi: 10. 1111/rda. 13844. [Epub ahead of print].

Yardimci, G. G., Ozadam, H., Sauria, M. E., Ursu, O., Yan, K.-K., Yang, T., et al. (2019). Measuring the reproducibility and quality of Hi-C data. Genome Biol. 20, 1–19. doi: 10. 1186/s13059-019-1658-7

Yin, Z., Ke, X., Li, Z., Chen, J., Gao, X., and Huang, L. (2017). Unconventional recombination in the mating type locus of heterothallic apple canker pathogen Valsa mali. G3 Genes Genomes Genet. 7, 1259–1265. doi: 10. 1534/g3. 116. 037853

Yu, A., Zhou, R., Xia, B., Dang, W., Yang, Z., and Chen, X. (2020). NAMPT maintains mitochondria content via NRF2-PPARα/AMPKα pathway to promote cell survival under oxidative stress. Cell. Signal 66: 109496. doi: 10. 1016/j. cellsig. 2019. 109496

Yuan, J., Liu, W., Karvar, S., Baker, S. S., He, W., Baker, R. D., et al. (2015). Potassium channel KCNJ15 is required for histamine-stimulated gastric acid secretion. Am. J. Physiol. Cell Physiol. 309, C264–C70. doi: 10. 1152/ajpcell. 00012. 2015

Zhan, T. T. Y., Lou, M., Liu, J., Wang, Y., and Zhao, X. (1994). The genetic mechanism of intersexuality in milk goats of Saanen breed of Xinong. Yi Chuan Xue Bao 21, 356–361.

Zhang, S., Cao, X., Li, Y., Wang, K., Yuan, M., and Lan, X. (2019). Detection of polled intersex syndrome (PIS) and its effect on phenotypic traits in goats. Anim. Biotechnol . 31, 561–565. doi: 10. 1080/10495398. 2019. 1625782

Zhang, Z., Chen, F., Li, S., Guo, H., Xi, H., Deng, J., et al. (2020). ERG the modulates Warburg effect and tumor progression in cervical cancer. Biochem. Biophys. Res. Commun. 522, 191–197. doi: 10. 1016/j. bbrc. 2019. 11. 079

Zheng, W., Nagaraju, G., Liu, Z., and Liu, K. (2012). Functional roles of the phosphatidylinositol 3-kinases (PI3Ks) signaling in the mammalian ovary. Mol. Cell. Endocrinol. 356, 24–30. doi: 10. 1016/j. mce. 2011. 05. 027

Zhou, X., Chen, Y., Mok, K. Y., Zhao, Q., Chen, K., Chen, Y., et al. (2018). Identification of genetic risk factors in the Chinese population implicates a role of immune system in Alzheimer’s disease pathogenesis. Proc. Natl. Acad. Sci. U. S. A . 115, 1697–1706. doi: 10. 1073/pnas. 1715554115