- Open Access
Identification of novel candidate genes for 46,XY disorders of sex development (DSD) using a C57BL/6J-Y POS mouse model
- Hayk Barseghyan1, 2,
- Aleisha Symon4,
- Mariam Zadikyan2,
- Miguel Almalvez1, 2,
- Eva E. Segura2,
- Ascia Eskin2,
- Matthew S. Bramble1, 2,
- Valerie A. Arboleda2,
- Ruth Baxter2,
- Stanley F. Nelson2,
- Emmanuèle C. Délot1, 2, 3,
- Vincent Harley†4 and
- Eric Vilain†1, 2, 3Email author
© The Author(s). 2018
- Received: 31 July 2017
- Accepted: 19 January 2018
- Published: 30 January 2018
Disorders of sex development (DSD) have an estimated frequency of 0.5% of live births encompassing a variety of urogenital anomalies ranging from mild hypospadias to a discrepancy between sex chromosomes and external genitalia. In order to identify the underlying genetic etiology, we had performed exome sequencing in a subset of DSD cases with 46,XY karyotype and were able to identify the causative genetic variant in 35% of cases. While the genetic etiology was not ascertained in more than half of the cases, a large number of variants of unknown clinical significance (VUS) were identified in those exomes.
To investigate the relevance of these VUS in regards to the patient’s phenotype, we utilized a mouse model in which the presence of a Y chromosome from the poschiavinus strain (Y POS ) on a C57BL/6J (B6) background results in XY undervirilization and sex reversal, a phenotype characteristic to a large subset of human 46,XY DSD cases. We assessed gene expression differences between B6-Y B6 and undervirilized B6-Y POS gonads at E11.5 and identified 515 differentially expressed genes (308 underexpressed and 207 overexpressed in B6-Y POS males).
We identified 15 novel candidate genes potentially involved in 46,XY DSD pathogenesis by filtering the list of human VUS-carrying genes provided by exome sequencing with the list of differentially expressed genes from B6-Y POS mouse model. Additionally, we identified that 7 of the 15 candidate genes were significantly underexpressed in the XY gonads of mice with suppressed Sox9 expression in Sertoli cells suggesting that some of the candidate genes may be downstream of a well-known sex determining gene, Sox9.
The use of a DSD-specific animal model improves variant interpretation by correlating human sequence variants with transcriptome variation.
- Disorders of sex development
- 46,XY DSD
- C57BL/6J mouse
- Gonadal dysgenesis
Human sex development is dictated by the inheritance of either an X or Y chromosome from the father to offspring. The male sex determination step starts with the expression of a Y-chromosome-encoded transcription factor SRY (sex-determining region Y) in the bipotential gonad, initiating a cascade of molecular and cellular events leading to testicular organogenesis . In the absence of the Y chromosome, female-specific pathways are initiated for proper ovarian development . Sex differentiation then occurs, mostly under the influence of testicular (e.g., testosterone, AMH) or ovarian (e.g., estradiol) hormones or transcription factors (e.g., COUP-TFII) that further differentiate the body into typical male or female structures, including both internal and external genitalia [3, 4]. Anomalies in hormonal exposure and/or gene mutations disrupting sex development pathways lead to disorders of sex development (DSD) [5–7], defined as “congenital conditions in which development of chromosomal, gonadal, or anatomic sex is atypical” . The umbrella term DSD encompasses conditions ranging from mild hypospadias (abnormal location of the meatus) to discrepancy between sex chromosomes and external genital phenotype (formerly known as sex reversal, either complete or with ambiguous genitalia). DSDs are estimated to affect up to 0.5% of the population . The birth of a child with a DSD may be highly stressful for families, bringing uncertainty in regard to the child’s future psychosexual development and clinical management [8, 10, 11].
At present, a specific molecular diagnosis is identified at variable rates in different DSD conditions, and gonadal dysgenesis cases are arguably the most difficult to diagnose. The majority (80–90%) of isolated 46,XX testicular DSD are explained by SRY translocations, but only a minority (~ 10%) of ovotesticular DSD in 46,XX individuals are . Copy number variants of the SOX9 and SOX3 gene regions are a well-established etiology but only explain a few cases . More recently, a single nucleotide variant in NR5A1 (nuclear receptor subfamily 5 group A member 1) gene resulting in p.Arg92Trp amino acid change has been associated with 46,XX testicular (and ovotesticular) DSD [13, 14]. The majority of cases of ovarian dysgenesis occur in individuals with an abnormal sex chromosome complement, most commonly 45,X (Turner syndrome), but the advent of next-generation sequencing has recently identified many autosomal genes implicated in determination and maintenance of the ovarian fate. They affect various processes, in particular DNA repair, replication, and stability, but explain a minority of cases [15, 16]. Among 46,XY DSD cases with gonadal dysgenesis, about 15% each are due to SRY, NR5A1, and MAP3K1 (mitogen-activated protein kinase kinase kinase 1), and rare cases have been attributed to mutations in other genes such as SOX9 (SRY-box9), NR0B1 (nuclear receptor subfamily 0 group B member 1), or FGFR2 (fibroblast growth factor receptor 2) [17, 18]. Nevertheless, collectively, the genetic etiology is still not identified in greater than 50% of DSD patients, suggesting the existence of a number of unknown sex-determining genes. We endeavored to identify novel candidate genes for 46,XY gonadal dysgenesis.
Next-generation sequencing has become instrumental in DSD diagnosis, including clinical exome sequencing and gene panels [17, 19–21] with high diagnostic rates reported for known DSD genes. In a cohort of 46,XY DSD patients, we established a diagnosis in approximately 1/3 of cases , similar to rates for other rare disorders [23, 24]. Another 15% of the exomes in the cohort contained variants of unknown significance (VUS) in known DSD genes that could not be validated as pathogenic but were reported to the referring clinicians to orient further endocrine or imaging testing toward a definitive diagnosis (the variants were termed as “actionable VUS”). Half of the cases from our cohort remained undiagnosed but contained hundreds of VUS that provide an opportunity for identification of novel etiologies for DSD. Here, we utilize an animal model of DSD with gonadal dysgenesis and undervirilization [25, 26] to identify a group of genes that were misexpressed during disrupted testis development. This list was cross-referenced with the list of VUS from 46,XY DSD patients to predict which VUS might be causative in cases where exome sequencing did not result in a definitive diagnosis. We show that the identified 15 novel candidate genes contain a VUS identified in 46,XY DSD cases and are expressed at the time of sex development in a sex-differential manner. In addition, we show that the expression of many of these genes in the developing male gonads is dependent on the known sex-determining gene Sox9.
Exome sequencing and analysis
DNA was isolated from peripheral blood using Gentra Puregene Blood Kit (Qiagen, USA) or saliva collected using ORAgene ORG-500 (DNAgenoteck, Canada). Sequencing libraries and exome capture was done for each sample following manufacturer’s protocols for SureSelect All Exon 50 Mb capture kit (Agilent Technologies) and Nextera Rapid Capture (Illumina, USA). Sequencing was performed on an Illumina HiSeq2500 as 50-100 bp paired-end run at the UCLA Clinical Genomics Center.
The sequence reads, FASTQ files, were aligned to the human reference genome (GRCh37/hg19 Feb. 2009 assembly) using BWA (Burrows-Wheeler Alignment tool)  and Novoalign (novocraft.com). The output BAM files were sorted and merged, and PCR duplicates were removed using Picard. INDEL (insertion and deletion) realignment and recalibration was performed using Genome Analysis Tool Kit (GATK) (broadinstitute.org). Both single-nucleotide variants (SNVs) and small INDELs were called within the Ensembl coding exonic intervals ± 2 bp using GATK’s Unified Genotyper, then recalibrated and filtered using GATK variant-quality score recalibration and variant filtration tools. All high-quality variants were annotated using SNP&Variation Suite and VarSeq—variant filtration and annotation software (Golden Helix, USA). All variants were filtered by a minor allele frequency (MAF) of < 1% and intersected with the DSD gene list to identify mutations in known DSD genes. The list is comprised of a primary gene list of well-annotated genes involved in sex determination and differentiation , as well as a secondary list of genes that are more loosely associated with sex development, e.g., their OMIM (Online Mendelian Inheritance of Man) description contains sex development keywords.
The variants identified by exome sequencing were classified into causative or likely causative variants following the recommendations of the American College of Medical Genetics and Genomics . All other variants with minor allele frequency below 1% were classified as variants of unknown significance (VUS). To assess previously unreported missense variants, we used two in silico algorithms SIFT  and PolyPhen  to predict the pathogenicity of a missense variant based on conservation of the amino acid across species, the physical characteristics of the altered amino acid, and the possible impact on protein structure and function. All variants with low quality scores were validated by Sanger sequencing .
Animal care and dissections
The C57BL/6J and C57BL/6J-Y POS animals were housed at the UCLA Animal Care Facility following the guidelines of the University of California, Los Angeles, Division of Laboratory Animal Medicine. All experiments were approved by the Institutional Animal Care and Research Committees of UCLA. Wild-type C57BL/6J males and females used for breeding were purchased from the Jackson Laboratory (Bar Harbor, USA), which is fully accredited by the American Association for Accreditation of Laboratory Animal Care.
We have previously identified a 1.5-Mb congenic region on chromosome 11 that confers 80% protection from B6-Y POS sex reversal in the heterozygous state (B6-110h-Y POS ) and complete protection in the homozygous state (B6-110H-Y POS ) . This protective region allows for continual maintenance of subfertile poschiavinus male mice as a breeding colony, with an option of generating unprotected B6-Y POS males by mating heterozygous B6-110h-Y POS males with wild-type (WT) B6 females. Overnight mating was performed using either the wild-type (WT) B6 or B6-110h-Y POS (protected from sex reversal) males and WT B6 females. Dissections were performed at E11.5; the gonads were separated from the mesonephros and placed in RNA stabilizing solution RNAlater (Ambion). DNA was extracted from the rest of the embryos for genotyping. Chromosomal sex was determined using a single primer pair for X-linked Smc-x gene (330 bp) and the Y-linked Smc-y gene (301 bp) (forward: 5′CCGCTGCCAAATTCTTTGG3′; reverse: 5′TGAAGCTTTTGGCTTTGAG3′). The presence of the Y POS chromosome was determined by a SNP between YB6 and Y POS Sry gene using the primer sets 5′TGAATGCATTTATGGTGTGGTC3′; 5′AGCTTTGCTGGTTTTTGGAGTA3′. Immomix Red (Bioline, UK). Presence or absence of the 110 h protective region in B6-Y POS males was checked by Sanger sequencing of two regions 11-10 and 11-11 containing SNP rs27019103 (5′AAAGTGTGCTTCCCAGGAGA3′; 5′CCTCTCCCTCAACCCCTAAG3′) and SNP rs28240850 (5′CCACAGCTGGAGGTAGGGTA3′; 5′CCTAAGATGCCATGGGAAGA3′) respectively .
Total RNA was isolated from combined embryonic gonadal tissue (50–70 gonads per group) using Qiagen RNeasy Kit (Qiagen) following manufacturer’s guidelines. RNA quality was assessed by Agilent 2100 Bioanalyzer (Agilent Technologies). All samples were required to have RNA integrity scores (RIN) greater than 8.
Experiments on AmhCre Sox9floxflox mice were carried out in strict adherence with the recommendations in the Australian code of practice for the care and use of animals for scientific purposes from the National Health and Medical Research Council.
E13.5 gonads were separated from the mesonephros and total RNA was isolated using the RNeasy Kit (Qiagen), as described in Rahmoun et al. . Embryos were genotyped, and the RNA from the six gonads was pooled into wild-type XY, XX, or XY AmhCre Sox9floxflox (Sox9 knockout). This was repeated in three biological replicates; protocols are detailed in Rahmoun et al. .
RNA sequencing and expression analysis
RNA from each sample was submitted to the UCLA Neuroscience Genomic Core (UNGC) for library preparation and sequencing. Library preparation was performed using TruSeq Stranded Total RNA kit (Illumina) with Poly-A selection following manufacturer’s guidelines. Sequencing was performed on HiSeq 2500 (Illumina) with 69 bp paired-end run on a rapid flowcell capable of generating 150 M reads per lane. Four samples were multiplexed and sequenced over two rapid lanes with each sample receiving approximately 75 million reads with > 85% map rate.
The generated sequencing reads were aligned to the mouse genome, version mm10 with STAR . Transcript abundance was assessed by Cufflinks (v2.1.1) , using a GTF file based on Ensembl mouse NCBI37. Differential expression analysis was based on fold change differences greater than 1.5 between the groups being compared. Differentially expressed genes were split into two categories: underexpressed and overexpressed in B6-Y POS males. Both categories were separately subjected to pathway enrichment analysis using Gene Ontology Consortium .
To analyze the RNA from AmhCre Sox9floxflox XY gonads and wild-type XY and XX gonads, libraries were generated using the NuGEN Mondrian Technology and SPIA amplification methodology, and the data was processed and aligned to the mouse genome (Ensembl version 38.77) as described by Rahmoun et al. . To eliminate composition biases, the trimmed mean of M values (TMM) method was used for normalization between the samples . The adjusted P value of 0.05 was used to assess which genes were differentially expressed between XY and AmhCre Sox9floxflox XY (Sox9 KO). Graphs of gene expression were made using GraphPad Prism.
Reverse transcription of RNA to cDNA was performed using Tetro cDNA Synthesis Kit (Bioline, UK) following manufacturer’s protocol. The primer sequences used are detailed in Additional file 1: Table S1. Primers were designed using autoprime software (autoprime.de) and spanned exon-exon junctions for optimal RNA quantification. cDNA was quantified using QuBit HS (Invitrogen) for double-stranded DNA, and a total of 3 ng of cDNA was used per sample for amplification. qPCR was carried out in duplicates using SensiFAST™ SYBR No-ROX Kit (Bioline, UK) by DNA Engine Opticon® 2 real-time PCR detection system (BioRad, USA). Reaction conditions were as follows: 95 °C for 10 min, then 40 cycles of 95 °C for 15 s, 60–64 °C (see Additional file 1: Table S1) for 10 s, and 72 °C for 15 s. Data was analyzed via Opticon Monitor Software (BioRad). Standard curves were generated from a mix of cDNA of all tested samples with five iterations of 1:4 dilutions. Average cycle threshold values (Ct) for each gene/sample were determined based on two replicates. Complementary DNA amounts were estimated based on Ct values and linear equation y = mx + b (where y is the Ct value, m is the slope, x is the cDNA amount, and b is the intercept).
Fbln2 expression in the embryonic gonads at E12.5 was assessed using immunohistochemistry following the experimental design of Wilhelm et al. , using the anti-Fbln2 rabbit polyclonal, sc-30176 (Santa Cruz Biotechnology) antibody. Topro (Invitrogen) was used to counterstain nuclei. All images were taken on a Zeiss LSM 510 Meta confocal microscope.
For the assessment of Sox9 and laminin expression in wild-type and AmhCre Sox9 floxflox gonads at E13.5, embryos were fixed overnight in 4% paraformaldehyde (PFA) at 4 °C, then washed three times in 1× PBS. The embryos were processed and embedded into paraffin, cut at 5 μm, and mounted onto slides. The slides were baked at 60 °C (30 min), deparaffinized using three washes of xylene, and hydrated using three washes of 100% ethanol, then distilled water and 1× PBS. Antigen retrieval was performed by microwaving slides (on high) in 10 mM sodium citrate (pH = 6.0) for 20 min. Sections were then blocked for 30 min with 5% normal donkey serum, and stained overnight at 4 °C with primary antibodies for anti-Sox9 rabbit polyclonal (1:400) and anti-Laminin rabbit polyclonal (1:100). Sections were washed three times in 1× PBS with 0.1% Tween20 (1× PBST) and incubated with the fluorescent-conjugated secondary antibodies, Donkey anti-rabbit AlexaFluor488 (Thermo Fisher, 1 μg/mL), for 1 h at room temp. Sections were washed three times in 1× PBST, then incubated in 0.1% Sudan Black in 70% EtOH for 5 min to quench background autofluorescence. Lastly, sections were washed three times in 1× PBST, counterstained using DAPI, then washed three times in 1× PBS and mounted using Dako Fluorescent Mounting Medium (Dako). Sections were imaged using fluorescence microscopy (Olympus Corp).
46,XY DSD cases with uninformative exome sequencing
Cohort of 46,XY DSD cases with uninformative exome sequencing
46,XY female, CGD
46,XY female, PGD
No uterus; Fallopian tubes present; short vagina; very low T and undetectable estradiol; gonads not found
46,XY female, GD
Amelia (missing limbs)
46,XY female, GD
Kidney disease; possible Denys-Drash syndrome
46,XY female, CGD
Normal uterus and Fallopian tubes; streak gonads
46,XY ambiguous genitalia
Partial fusion of labioscrotal folds; small phallus; penoscrotal hypospadias
46,XY ambiguous genitalia
Developmental delay; agenesis of corpus callosum
46,XY ambiguous genitalia
Adrenal hypoplasia congenita; dysmorphic features
46,XY ambiguous genitalia
Microcephaly; intestinal dysmotility; optic nerve hypoplasia
46,XY male, micropenis/cryptochidism
Severe growth and developmental retardation; testes not seen by ultrasound
46,XY male, hypospadias
Large clitoris; no uterus or vaginal opening; inguinal testes
46,XY ambiguous genitalia, CGD
Abdominal gonads with no oocytes; no seminiferous tubules; no clitoromegaly; posterior fusion of labia; urogenital sinus
Inguinal testes w/ immature seminiferous tubules; no uterus or Fallopian tubes; deafness; impaired cognition
46,XY ambiguous genitalia
Undescended testes; bifid scrotum; hypospadias
46,XY ambiguous genitalia
Bilateral descended testes; midshaft hypospadias; chordee
46,XY male, micropenis
No uterus or ovaries per ultrasound; ambiguous genitalia; undervirilization
Complete androgen insensitivity syndrome
46,XY male, hypospadias
46,XY female, GD
46,XY male, anorchia
Congenital bilateral anorchia; fully formed scrotum; definite penis (mildly shortened); no hypospadias; responsive to testosterone
46,XY male, hypospadias/cryptorchidism
Azoospermia; high T levels
Multiple congenital anomalies; no uterus; abdominal gonads—testes
46,XY male, microphallus
46,XY male, micropenis
46,XY male, hypospadias
46,XY male, hypospadias
Chordee; bifid scrotum; cryptorchidism
Growth delay; short stature
C57BL/6J-Y poschiavinus mice as a model for 46,XY gonadal dysgenesis
In cases where no pathogenic variant was found by exome sequencing, we identified many VUS outside of the DSD clinical gene list. To investigate the relevance of these VUS in regard to patients’ phenotype, we utilized a powerful mouse model for studying undervirilization in human 46,XY individuals. In this model, the presence of a Y chromosome originating from a M. domesticus poschiavinus strain (Y POS ) on a C57BL/6J (B6) background (B6-Y POS ), an inbred laboratory strain that normally carries a M. musculus Y chromosome, results in disrupted testicular development and female genital phenotype .
In the embryonic mouse gonad, Sry is normally expressed in a dynamic wave (central to distal) between E10.5 and E12.5 in the XY genital ridge, with peak Sry expression occurring in normal XY B6 genital ridges at ~E11.5, i.e., at the 16–18 tail somite stage of development, which is followed by the upregulation of Sox9 [39, 40]. In contrast, expression of the Sry POS gene peaks 10 to 14 h later in the genital ridges of B6-Y POS fetuses . We hypothesized that abnormal gonadal expression of specific genes in B6-Y POS males, after the surge of Sry during gonadal development, would correlate with the genes in which VUS were identified in 46,XY DSD patients by exome sequencing.
Gene expression differences between B6-Y B6 and of B6-Y POS males
Since all of the 46,XY DSD patients in the cohort carried a functional SRY gene, it was important to perform gene expression analysis in the animal model after the peak of Sry expression for optimal comparability. To achieve this, gonadal tissue from WT B6-Y B6 and undervirilized B6-Y POS males at embryonic day E11.5, specifically at 21 tail somites (a time point when the surge of Sry gene was complete in both B6-Y B6 and B6-Y POS males), was collected to perform RNA sequencing for assessment of differential gene expression.
Using this method, we identified 515 genes that were differentially expressed between B6-Y B6 and B6-Y POS males with a fold change greater than 1.5. Out of these 515 genes, 308 were underexpressed and 207 were overexpressed in B6-Y POS males (Fig. 1b; Additional file 2: Table S2). To validate the integrity of the method and make sure that the correct tissue was dissected at the correct embryonic stages, we first looked at the expression levels of two important genes involved in testicular development, Sry and Sox9. High Sry and low Sox9 expression levels in B6-Y POS males indicated the correct timing of embryonic development (Fig. 1c), expression of which coincided with a previous publication . Second, we verified which genes in our DSD gene list used for exome variant filtering were present in the B6-Y B6 /B6-Y POS differentially expressed gene list. The comparison of the two lists revealed that 21 genes were in common: 15 underexpressed and 6 overexpressed (Additional file 2: Table S2). In our previous cohort , out of these 21 genes, 3 (HSD17B3 (hydroxysteroid 17-beta dehydrogenase 3), STAR (steroidogenic acute regulatory protein), FGFR2) contained a pathogenic variant identified by exome sequencing, explaining a total of 5 cases, and 2 (DHH (desert hedgehog), MAMLD1 (mastermind-like domain-containing 1)) contained a variant that was reported to the clinician to orient further endocrine or imaging testing toward a definitive diagnosis. Cumulatively, these findings indicate that the gene expression analyses were carried out in a correct tissue type, at the correct developmental time point, and that the differentially expressed genes between B6-Y B6 and B6-Y POS males may potentially be important in sex development.
Filtering of VUS in 46,XY DSD cases using the B6-Y POS gene list
On average, exome sequencing identifies ~ 21,000 variants per single case . Since DSDs are rare conditions, all variants identified in exome with a minor allele frequency (MAF) of more than 1% in the population were excluded. The variants remaining after the MAF cutoff were classified as variants of unknown significance. The number of genes with a VUS in each case ranged from 30 to 1100 with an average of approximately 730 genes per case. The gene list generated via expression studies in B6-Y B6 /B6-Y POS mice, consisting of 515 genes, was used to filter the list of VUS-containing human genes identified by exome sequencing. The comparison of two lists identified 305 (189 underexpressed and 116 overexpressed in B6-Y POS ) genes that were both differentially expressed in B6-Y POS males and contained a VUS in the 46,XY DSD cohort with an uninformative exome (Fig. 1b).
All these genes are known to be expressed in the developing gonad at the time of sex determination (e.g., the method used to identify these genes intrinsically already ensures that all those genes are expressed in the relevant tissue (gonad) at the relevant developmental time). In order to increase the probability of identifying relevant candidate genes involved in 46,XY DSD pathogenesis, we further queried if the differentially expressed genes from the mouse model (all 515 genes) were involved in any known biological processes. Gene Ontology Consortium (GOC)  enrichment analysis confirmed that genes underexpressed in B6-Y POS males were indeed enriched in biological processes known to control multicellular organism and anatomical structure development, including male reproductive development (Additional file 3: Table S3). Understanding the relevance of the genes that were overexpressed in B6-Y POS males was less straightforward. These genes were enriched in only two biological processes: response to extracellular stimulus and epithelial cell differentiation. Both of these categories had a high P value indicating that many genes in the overexpressed category are not associated with any known biological processes at this time. In addition, all of the pathogenic variants identified in our previous 46,XY DSD cohort  were in the underexpressed category of genes, indicating that they need to be expressed at higher levels in the developing gonad for proper testicular formation.
Based on these findings, we focused on variants identified in genes underexpressed in B6-Y POS males whose higher expression in WT males correlated with normal male sex development. To choose a fold change cutoff, we looked at fold change differences in expression between B6-Y B6 and B6-Y POS males for genes present in our clinical primary gene list. We found that the majority of the genes have an expression that is 2-fold higher in WT males compared in B6-Y POS males (Fig. 1d). To make the analysis more stringent and improve the confidence of identifying true candidate genes involved in male sex development, we therefore increased the fold change cutoff in expression between the B6-Y B6 and B6-Y POS males from 1.5 to 2. This change decreased the number of underexpressed genes from 189 to 53. Additional filtering was performed based on variant frequency (variants with MAF close to or below 0.1%), amino acid conservation (variants in highly conserved residues across multiple species were given preference), number of variants contained in a gene across the cohort, in silico predictions for pathogenicity (preference was given to the variants predicted to be deleterious or damaging), availability of literature (some weight was given to genes with known functions), and gonadal cell-specific expressivity using GenitoUrinary Developmental Molecular Anatomy Project (GUDMAP) data  (preference was given to genes expressed in male-typical cells).
List of VUS in candidate genes found in the cohort of 32 46,XY DSD patients
DSD case ID
MAF gnomAD (%)
Expression of the novel candidate genes is Sox9-dependent
The use of the undervirilized B6-Y POS mice as a model for 46,XY DSD in humans provides valuable screening information toward the identification of novel genes involved in male sex development, mutations in which could lead to anomalies in gonadal development in 46,XY patients with DSD. All of the identified candidate genes are expressed in the developing mouse gonad at the relevant time for sex determination and, as we have shown (Fig. 3a, b), the expression of many of these genes may be Sox9-dependent. However, when studying complex disorders such as DSD, it is important to note that even though the mouse models used here are extremely beneficial for identification of the underlying genetic cause in humans, they still do not provide the full spectrum of gene expression/interactions that occur during human sex development.
Mutations in the novel candidate genes identified via the Y POS mouse model are likely to be causative. For example, one of the candidate genes Adamts16 (A disintegrin and metalloproteinase with thrombospondin type 1 motif, 16) has been shown to be co-expressed with the known DSD gene Wt1 (Wilms tumor 1) in embryonic gonads, adult testes, and spermatids . Moreover, targeted disruption of Adamts16 in rats results in cryptorchidism and sterility . In our 46,XY DSD cohort, we identified three heterozygous variants in this gene (Table 2). Patients RDSD013 and RDSD002, both 46,XY women with complete gonadal dysgenesis, had a missense variant leading to amino acid change at positions p.Val734Ile and p.Arg100Trp respectively. These changes were located in the propeptide or cysteine-rich domain of the ADAMTS16 protein and may prevent expression or proper folding of the protein. The third missense variant (p.Phe469Val) in patient RDSD022 (46,XY, with ambiguous genitalia) was located in the peptidase domain of the protein and predicted damaging by in silico tools suggesting a possible impairment of the enzymatic function of ADAMTS16.
We identified a single variant, predicted to be damaging by in silico tools, in the SPRY4 (sprouty RTK signaling antagonist 4) gene in a 46,XY male patient (CDSD039) with hypogonadism. SPRY4 variants have been found in a cohort of patients presenting with hypogonadotropic hypogonadism with or without anosmia (HH17, OMIM #615266) . These genes are believed to be functioning in an oligogenic model, with variants in several genes possibly needed for phenotypic expression. Variants in SPRY4 have been found in association with variants in FGFR1 (fibroblast growth factor receptor 1) (HH2, OMIM #147950) and DUSP6 (dual specificity phosphatase 6) (HH19, OMIM #615269), the two other FGF signaling pathway components. An FGFR2 missense mutation was reported in a 46,XY female DSD patient, for which a corresponding mouse model showed partial sex reversal with reduced Spry4 (2-fold) and Dusp6 expression (> 2-fold) . We did not identify FGFR1 or DUSP6 variants in the exome of patient CDSD039 (which would have been diagnostic for this patient). However, DUSP6 is present in the differentially expressed gene list (underexpressed in B6-Y POS with a fold change of 1.7) and another gene coding for a dual-specificity phosphatase, DUSP15 (dual specificity phosphatase 15), is in our final candidate gene list, with underexpression in B6-Y POS (fold change > 2) and contains a VUS in one patient.
Exome sequencing provides high-throughput genetic diagnostic capability that has become the core of modern clinical genetics. However, many variants identified by whole exome sequencing are uninterpretable clinically. The C57/BL6J-Y POS model narrows the interpretive gap by correlating human sequence variants with transcriptome variation. This approach allowed the identification of 15 novel candidate genes for human 46,XY DSD.
We would like to thank all the patients and families who contributed samples to this project and consented to have their data shared. We would also like to thank the clinical and research teams of the DSD-Translational Research Network, founded by Eric Vilain (UCLA) and David E. Sandberg (University of Michigan, Ann Arbor) with the collaborations of DSD patient advocate organization Accord Alliance.
This work was supported by the Doris Duke Foundation (to E.V.), a National Institutes of Health T032 Training Grant (5T32HD007228) (to H.B.), Eunice Kennedy Shriver National Institute of Child Health and Human Development (NICHD) Grant (RO1HD06138), Disorders of Sex Development-Translational Research Network (to EV), and National Health and Medical Research Council Program Grants 334314 and 546517 (to V.H.) and Fellowships 441102 and 1020034 (to V.H.), and Australian Government Research Training Program Scholarship (to A.S.). This work was also supported by the Victorian Government’s Operational Infrastructure Support Program.
Availability of data and materials
Data generated or analyzed during this study are partially included in this published article. The datasets generated and/or analyzed during the current study are available from the corresponding author on a reasonable request.
HB designed the study; acquired, analyzed, and interpreted the data; and prepared the manuscript. AS, MZ, MA, ES, AS, MSB, VAA, and RB acquired, analyzed, and interpreted the data. AS, ED, VH, and EV aided in the experimental design, data analysis, and manuscript preparation. SFN, ED, VH, and EV provided senior oversight for the design, analysis, and interpretation. All authors indicated the approval for publication.
Ethics approval and consent to participate
Research involving human subjects, human material, and human data have been performed in accordance with protocols (IRB# 11-001491; IRB# 11-001775) approved by the UCLA Institutional Review Board.
Experimental studies involving animals, animal tissues, and data have been performed in accordance with approved protocol (ARC# 2000-088-51A) by UCLA Chancellor’s Animal Research Committee as well as Animal Ethics Committee of Southern Health (Clayton, Australia, ethics number MMCB/2009/30).
Consent for publication
Participants have consented by approved IRB protocol to share their de-identified data.
The authors declare that they have no competing interests.
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
- Koopman P, Sinclair A, Lovell-Badge R. Of sex and determination: marking 25 years of Randy, the sex-reversed mouse. Development. 2016;143(10):1633–7.View ArticlePubMedGoogle Scholar
- Eggers S, Ohnesorg T, Sinclair A. Genetic regulation of mammalian gonad development. Nat Rev Endocrinol. 2014;10(11):673–83.View ArticlePubMedGoogle Scholar
- Ahmed SF, Hughes IA. The genetics of male undermasculinization. Clin Endocrinol. 2002;56(1):1–18.View ArticleGoogle Scholar
- Zhao F, Franco HL, Rodriguez KF, Brown PR, Tsai MJ, Tsai SY, et al. Elimination of the male reproductive tract in the female embryo is promoted by COUP-TFII in mice. Science. 2017;357(6352):717–20.View ArticlePubMedPubMed CentralGoogle Scholar
- Arboleda VA, Fleming AA, Vilain E. Disorders of sex development. In: Weiss RE, Refetoff S, editors. Genetic diagnosis of endocrine disorders. London: Academic Press; 2010. p. 227–43.Google Scholar
- Délot E, Vilain E. Disorders of sex development. In: Strauss JF, Barbieri RL, Gargiulo AR, editors. Yen & Jaffe’s reproductive endocrinology. Philadelphia: Elsevier; 2019.Google Scholar
- Ono M, Harley VR. Disorders of sex development: new genes, new concepts. Nat Rev Endocrinol. 2013;9(2):79–91.View ArticlePubMedGoogle Scholar
- Lee PA, Houk CP, Ahmed SF, Hughes IA. Consensus statement on management of intersex disorders. International Consensus Conference on Intersex. Pediatrics. 2006;118(2):e488–500.View ArticlePubMedGoogle Scholar
- Lee PA, Nordenstrom A, Houk CP, Ahmed SF, Auchus R, Baratz A, et al. Global disorders of sex development update since 2006: perceptions, approach and care. Horm Res Paediatr. 2016;85(3):158–80.View ArticlePubMedGoogle Scholar
- Sandberg DE, Gardner M, Cohen-Kettenis PT. Psychological aspects of the treatment of patients with disorders of sex development. Semin Reprod Med. 2012;30(5):443–52.View ArticlePubMedPubMed CentralGoogle Scholar
- Warne GL. Long-term outcome of disorders of sex development. Sex Dev. 2008;2(4–5):268–77.View ArticlePubMedGoogle Scholar
- Delot EC, Vilain EJ. Nonsyndromic 46,XX testicular disorders of sex development. In: Pagon RA, et al., editors. GeneReviews(R). Seattle; 1993.Google Scholar
- Baetens D, Stoop H, Peelman F, Todeschini AL, Rosseel T, Coppieters F, et al. NR5A1 is a novel disease gene for 46,XX testicular and ovotesticular disorders of sex development. Genet Med. 2017;19(4):367–76.View ArticlePubMedGoogle Scholar
- Bashamboo A, Donohoue PA, Vilain E, Rojo S, Calvel P, Seneviratne SN, et al. A recurrent p.Arg92Trp variant in steroidogenic factor-1 (NR5A1) can act as a molecular switch in human sex development. Hum Mol Genet. 2016;25(23):5286.PubMedGoogle Scholar
- Chapman C, Cree L, Shelling AN. The genetics of premature ovarian failure: current perspectives. Int J Womens Health. 2015;7:799–810.PubMedPubMed CentralGoogle Scholar
- Caburet S, Arboleda VA, Llano E, Overbeek PA, Barbero JL, Oka K, et al. Mutant cohesin in premature ovarian failure. N Engl J Med. 2014;370(10):943–9.View ArticlePubMedPubMed CentralGoogle Scholar
- Barseghyan H, Delot E, Vilain E. New genomic technologies: an aid for diagnosis of disorders of sex development. Horm Metab Res. 2015;47(5):312–20.View ArticlePubMedGoogle Scholar
- Granados A, Alaniz VI, Mohnach L, Barseghyan H, Vilain E, Ostrer H, et al. MAP3K1-related gonadal dysgenesis: six new cases and review of the literature. Am J Med Genet C Semin Med Genet. 2017;175(2):253–9.View ArticlePubMedGoogle Scholar
- Delot EC, Papp JC, DSD-TRN Genetics Workgroup, Sandberg DE, Vilain E. Genetics of disorders of sex development: the DSD-TRN experience. Endocrinol Metab Clin N Am. 2017;46(2):519–37.View ArticleGoogle Scholar
- Eggers S, Sadedin S, van den Bergen JA, Robevska G, Ohnesorg T, Hewitt J, et al. Disorders of sex development: insights from targeted gene sequencing of a large international patient cohort. Genome Biol. 2016;17(1):243.View ArticlePubMedPubMed CentralGoogle Scholar
- Kim JH, Kang E, Heo SH, Kim GH, Jang JH, Cho EH, et al. Diagnostic yield of targeted gene panel sequencing to identify the genetic etiology of disorders of sex development. Mol Cell Endocrinol. 2017;444:19–25.View ArticlePubMedGoogle Scholar
- Baxter RM, Arboleda VA, Lee H, Barseghyan H, Adam MP, Fechner PY, et al. Exome sequencing for the diagnosis of 46,XY disorders of sex development. J Clin Endocrinol Metab. 2015;100(2):E333–44.View ArticlePubMedGoogle Scholar
- Lee H, Deignan JL, Dorrani N, Strom SP, Kantarci S, Quintero-Rivera F, et al. Clinical exome sequencing for genetic identification of rare Mendelian disorders. JAMA. 2014;312(18):1880–7.View ArticlePubMedPubMed CentralGoogle Scholar
- Yang Y, Muzny DM, Reid JG, Bainbridge MN, Willis A, Ward PA, et al. Clinical whole-exome sequencing for the diagnosis of mendelian disorders. N Engl J Med. 2013;369(16):1502–11.View ArticlePubMedPubMed CentralGoogle Scholar
- Arboleda VA, Fleming A, Barseghyan H, Delot E, Sinsheimer JS, Vilain E. Regulation of sex determination in mice by a non-coding genomic region. Genetics. 2014;197(3):885–97.View ArticlePubMedPubMed CentralGoogle Scholar
- Umemura Y, Miyamoto R, Hashimoto R, Kinoshita K, Omotehara T, Nagahara D, et al. Ontogenic and morphological study of gonadal formation in genetically-modified sex reversal XY(POS) mice. J Vet Med Sci. 2016;77(12):1587–98.View ArticlePubMedGoogle Scholar
- Li H, Durbin R. Fast and accurate long-read alignment with burrows-wheeler transform. Bioinformatics. 2010;26(5):589–95.View ArticlePubMedPubMed CentralGoogle Scholar
- Rehm HL, Bale SJ, Bayrak-Toydemir P, Berg JS, Brown KK, Deignan JL, et al. ACMG clinical laboratory standards for next-generation sequencing. Genet Med. 2013;15(9):733–47.View ArticlePubMedPubMed CentralGoogle Scholar
- Kumar P, Henikoff S, Ng PC. Predicting the effects of coding non-synonymous variants on protein function using the SIFT algorithm. Nat Protoc. 2009;4(7):1073–81.View ArticlePubMedGoogle Scholar
- Adzhubei IA, Schmidt S, Peshkin L, Ramensky VE, Gerasimova A, Bork P, et al. A method and server for predicting damaging missense mutations. Nat Methods. 2010;7(4):248–9.View ArticlePubMedPubMed CentralGoogle Scholar
- Strom SP, Lee H, Das K, Vilain E, Nelson SF, Grody WW, et al. Assessing the necessity of confirmatory testing for exome-sequencing results in a clinical molecular diagnostic laboratory. Genet Med. 2014;16(7):510–5.View ArticlePubMedPubMed CentralGoogle Scholar
- Rahmoun M, Lavery R, Laurent-Chaballier S, Bellora N, Philip GK, Rossitto M, et al. In mammalian foetal testes, SOX9 regulates expression of its target genes by binding to genomic regions with conserved signatures. Nucleic Acids Res. 2017;45(12):7191–211.View ArticlePubMedPubMed CentralGoogle Scholar
- Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29(1):15–21.View ArticlePubMedGoogle Scholar
- Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–5.View ArticlePubMedPubMed CentralGoogle Scholar
- The Gene Ontology, C. Expansion of the Gene Ontology knowledgebase and resources. Nucleic Acids Res. 2017;45(D1):D331–8.View ArticleGoogle Scholar
- Robinson MD, Oshlack A. A scaling normalization method for differential expression analysis of RNA-seq data. Genome Biol. 2010;11(3):R25.View ArticlePubMedPubMed CentralGoogle Scholar
- Wilhelm D, Washburn LL, Truong V, Fellous M, Eicher EM, Koopman P. Antagonism of the testis- and ovary-determining pathways during ovotestis development in mice. Mech Dev. 2009;126(5–6):324–36.View ArticlePubMedPubMed CentralGoogle Scholar
- Eicher EM, Washburn LL, Whitney JB 3rd, Morrow KE. Mus poschiavinus Y chromosome in the C57BL/6J murine genome causes sex reversal. Science. 1982;217(4559):535–7.View ArticlePubMedGoogle Scholar
- Bullejos M, Koopman P. Spatially dynamic expression of Sry in mouse genital ridges. Dev Dyn. 2001;221(2):201–5.View ArticlePubMedGoogle Scholar
- Kobayashi A, Chang H, Chaboissier MC, Schedl A, Behringer RR. Sox9 in testis determination. Ann N Y Acad Sci. 2005;1061:9–17.View ArticlePubMedGoogle Scholar
- Bullejos M, Koopman P. Delayed Sry and Sox9 expression in developing mouse gonads underlies B6-Y(DOM) sex reversal. Dev Biol. 2005;278(2):473–81.View ArticlePubMedGoogle Scholar
- Jameson SA, Natarajan A, Cool J, DeFalco T, Maatouk DM, Mork L, et al. Temporal transcriptional profiling of somatic and germ cells reveals biased lineage priming of sexual fate in the fetal mouse gonad. PLoS Genet. 2012;8(3):e1002575.View ArticlePubMedPubMed CentralGoogle Scholar
- Harding SD, Armit C, Armstrong J, Brennan J, Cheng Y, Haggarty B, et al. The GUDMAP database—an online resource for genitourinary research. Development. 2011;138(13):2845–53.View ArticlePubMedPubMed CentralGoogle Scholar
- Barrionuevo F, Georg I, Scherthan H, Lecureuil C, Guillou F, Wegner M, et al. Testis cord differentiation after the sex determination stage is independent of Sox9 but fails in the combined absence of Sox9 and Sox8. Dev Biol. 2009;327(2):301–12.View ArticlePubMedGoogle Scholar
- Barrionuevo F, Bagheri-Fam S, Klattig J, Kist R, Taketo MM, Englert C, et al. Homozygous inactivation of Sox9 causes complete XY sex reversal in mice. Biol Reprod. 2006;74(1):195–201.View ArticlePubMedGoogle Scholar
- Lavery R, Lardenois A, Ranc-Jianmotamedi F, Pauper E, Gregoire EP, Vigier C, et al. XY Sox9 embryonic loss-of-function mouse mutants show complete sex reversal and produce partially fertile XY oocytes. Dev Biol. 2011;354(1):111–22.View ArticlePubMedGoogle Scholar
- Bi W, Huang W, Whitworth DJ, Deng JM, Zhang Z, Behringer RR, et al. Haploinsufficiency of Sox9 results in defective cartilage primordia and premature skeletal mineralization. Proc Natl Acad Sci U S A. 2001;98(12):6698–703.View ArticlePubMedPubMed CentralGoogle Scholar
- Behringer RR, Finegold MJ, Cate RL. Müllerian-inhibiting substance function during mammalian sexual development. Cell. 1994;79(3):415–25.View ArticlePubMedGoogle Scholar
- Moniot B, Declosmenil F, Barrionuevo F, Scherer G, Aritake K, Malki S, et al. The PGD2 pathway, independently of FGF9, amplifies SOX9 activity in Sertoli cells during male sexual differentiation. Development. 2009;136(11):1813–21.View ArticlePubMedPubMed CentralGoogle Scholar
- Kashimada K, Svingen T, Feng CW, Pelosi E, Bagheri-Fam S, Harley VR, et al. Antagonistic regulation of Cyp26b1 by transcription factors SOX9/SF1 and FOXL2 during gonadal development in mice. FASEB J. 2011;25(10):3561–9.View ArticlePubMedPubMed CentralGoogle Scholar
- Abdul-Majeed S, Mell B, Nauli SM, Joe B. Cryptorchidism and infertility in rats with targeted disruption of the Adamts16 locus. PLoS One. 2014;9(7):e100967.View ArticlePubMedPubMed CentralGoogle Scholar
- Bouma GJ, Hudson QJ, Washburn LL, Eicher EM. New candidate genes identified for controlling mouse gonadal sex determination and the early stages of granulosa and Sertoli cell differentiation. Biol Reprod. 2010;82(2):380–9.View ArticlePubMedGoogle Scholar
- Miraoui H, Dwyer AA, Sykiotis GP, Plummer L, Chung W, Feng B, et al. Mutations in FGF17, IL17RD, DUSP6, SPRY4, and FLRT3 are identified in individuals with congenital hypogonadotropic hypogonadism. Am J Hum Genet. 2013;92(5):725–43.View ArticlePubMedPubMed CentralGoogle Scholar
- Bagheri-Fam S, Ono M, Li L, Zhao L, Ryan J, Lai R, et al. FGFR2 mutation in 46,XY sex reversal with craniosynostosis. Hum Mol Genet. 2015;24(23):6699–710.View ArticlePubMedPubMed CentralGoogle Scholar
- Hennekam RC, Allanson JE, Biesecker LG, Carey JC, Opitz JM, Vilain E. Elements of morphology: standard terminology for the external genitalia. Am J Med Genet A. 2013;161A(6):1238–63.View ArticlePubMedGoogle Scholar