GnRH Deficient Patients With Congenital Hypogonadotropic Hypogonadism: Novel Genetic Findings in ANOS1, RNF216, WDR11, FGFR1, CHD7, and POLR3A Genes in a Case Series and Review of the Literature

Background: Congenital hypogonadotropic hypogonadism (CHH) is a rare genetic disease caused by Gonadotropin-Releasing Hormone (GnRH) deficiency. So far a limited number of variants in several genes have been associated with the pathogenesis of the disease. In this original research and review manuscript the retrospective analysis of known variants in ANOS1 (KAL1), RNF216, WDR11, FGFR1, CHD7, and POLR3A genes is described, along with novel variants identified in patients with CHH by the present study. Methods: Seven GnRH deficient unrelated Cypriot patients underwent whole exome sequencing (WES) by Next Generation Sequencing (NGS). The identified novel variants were initially examined by in silico computational algorithms and structural analysis of their predicted pathogenicity at the protein level was confirmed. Results: In four non-related GnRH males, a novel X-linked pathogenic variant in ANOS1 gene, two novel autosomal dominant (AD) probably pathogenic variants in WDR11 and FGFR1 genes and one rare AD probably pathogenic variant in CHD7 gene were identified. A rare autosomal recessive (AR) variant in the SRA1 gene was identified in homozygosity in a female patient, whilst two other male patients were also, respectively, found to carry novel or previously reported rare pathogenic variants in more than one genes; FGFR1/POLR3A and SRA1/RNF216. Conclusion: This report embraces the description of novel and previously reported rare pathogenic variants in a series of genes known to be implicated in the biological development of CHH. Notably, patients with CHH can harbor pathogenic rare variants in more than one gene which raises the hypothesis of locus-locus interactions providing evidence for digenic inheritance. The identification of such aberrations by NGS can be very informative for the management and future planning of these patients.


INTRODUCTION
Congenital hypogonadotropic hypogonadism (CHH) is a rare disorder that is mainly caused by gonadotropin releasing hormone (GnRH) deficiency and characterized by delayed sexual development and infertility in both males and females (1)(2)(3)(4)(5)(6). The pulsatile secretion of the decapeptide GnRH from the hypothalamus into the hypophyseal-portal vessels exerts control in the synthesis and release of luteinizing hormone (LH) and follicle stimulating hormone (FSH) in the anterior pituitary gland (7,8). In CHH, the GnRH secretion and/or action is impaired and as a consequence patients with the disorder exhibit low levels of gonadotropins, low sex steroids, absent, incomplete or delayed puberty and subsequently hypogonadotropic hypogonadism (HH) (1,9,10). The prevalence of CHH is estimated to be 1:8,000 male and 1:40,000 female live births with slightly fewer than 50% of cases suffering from hyposmia or anosmia (10)(11)(12). CHH is divided into two subtypes, which include of congenital normosmic isolated hypogonadotropic hypogonadism (HH) and anosmic HH or Kallmann syndrome (KS) (13,14). Current research regarding the pathophysiology of CHH provides evidence that genetic abnormalities play a key role in the development of the disease and is estimated that a genetic cause is apparent in almost 50% of CHH cases (1,6). Up-to-date there have been reported more than 60 putative loci for CHH, 17 of which have been linked with KS (1,6,(13)(14)(15)(16).
Over the last few years with the use of the high throughput next generation sequencing (NGS) the number of genes shown to be responsible for causing CHH/KS has radically increased (15,17). Therefore, the purpose of the present study was to determine the genetic involvement in a series of clinically diagnosed with CHH/KS Cypriot patients.

MATERIALS AND METHODS Patients
A total of seven (six males and one female) unrelated Cypriot patients with CHH/KS were included in the present study and underwent whole exome sequencing by NGS. Clinical criteria included the absence or incomplete development of secondary sexual characteristics after the age of 16 years in females and 18 years in males. The biochemical criteria included low levels of basal and GnRH stimulated gonadotrophins (LH, FSH) as well as low levels of sex steroids (testosterone in males and estradiol in females). MRI scans were performed for all patients, with the exception of patients 2 and 6. Five male patients and one female exhibited isolated hypogonadotropic hypogonadism. Only one patient, a 72-year-old male, had KS with anosmia. Written, informed consent was obtained from all seven adult individuals that participate in the study for the publication of any potentially identifiable images or data included in this article. The study was approved by the Cyprus National Bioethics Committee and all methods were performed in accordance with the relevant guidelines and regulations.

Genetic Analysis
Genomic DNA was isolated from peripheral blood using the Gentra Puregene Kit (Qiagen, Valencia, CA, USA) according to the manufacturer's instructions. The DNA concentration and purity was measured using the Nanodrop ND-1000 spectrophotometer (NanoDrop Technologies, Wilmington, DE, USA). Prior to library preparation for whole exome sequencing (WES) genomic DNA was quantified using the Qubit dsDNA BR Assay Kit (Invitrogen, Life Technologies, Eugene, OR, USA) on a Qubit R 2.0 Fluorometer (Invitrogen, Life Technologies, Eugene, OR, USA). WES was performed by using the TruSeq Exome Kit (Illumina Inc., San Diego, CA, USA) with pairedend 150 bp reads. NGS was performed using the NextSeq 500/550 High Output Kit v2.5 (150 Cycles) on an NextSeq500 system (Illumina Inc., San Diego, CA, USA). The FastQC quality control tool (http://www.bioinformatics.babraham.ac.uk/ projects/fastqc/) was used to evaluate the quality of the WES procedure. The mean target coverage for the exome was 70.88X. Specifically, 10X coverage was achieved for 99.22% of the nucleotides, 20X coverage for 87.68% of the nucleotides and 30X coverage for 79.35% of the nucleotides, indicating that the WES reaction was of sufficiently high quality for subsequent analysis.

Variant Analysis
The fastqc data obtained by WES were processed using an in-house bioinformatics pipeline. Briefly, all variants were inputted into the VarApp Browser and filtered. VarApp is a graphical user interface, which supports GEMINI (18). Variants in selected genes with biological involvement in the GnRH neuronal system and CHH (Supplementary Table 1) were further analyzed using the Qualimap v2.2.1 tool (19) to calculate the target coverage. Mean target coverage was >20X for 93.2% of the selected genes and >30X for 89% of the selected genes (Supplementary Table 1). Variants in these genes were additionally filtered using the VarApp Browser for minor allele frequencies of <1% in public databases such as 1000 genomes, ExAC browser and Exome Sequencing Project (ESP). Moreover, variants were filtered and selected according to their impact such as frameshift, splice acceptor, splice donor, start lost, stop gained, stop lost, inframe deletion, inframe insertion, missense, protein altering and splice region. In addition, variants were filtered by the VarApp Browser for their pathogenicity by two in silico tools, SIFT and Polyphen2. Population-specific data from an in-house WES library composed of 51 randomly selected samples of Cypriot origin were used to evaluate the potential disease-causing variants. All variants identified were confirmed by Sanger sequencing. When genetic material of relatives was obtainable, familial segregation was performed (Table 2, Supplementary Figure 5). For the cases where two potentially pathogenic variants were identified in an individual, we employed the ORVAL platform for predicting pathogenicity due to digenicity (20). Finally, the variants were categorized for their pathogenicity using the standards and guidelines of the American College of Medical Genetics and Genomics and the Association for Molecular Pathology (21).

In silico Analysis of the Single Nucleotide Variants
In silico prediction on protein function of the pathogenicity effect of the different amino acid substitutions identified by NGS and confirmed with Sanger analysis was performed by the PredictSNP tool using the default settings (22). The PredictSNP tool evaluates the pathogenicity of a variant by using seven different protein functionality prediction tools: MAPP, PhD-SNP, PolyPhen-1, PolyPhen-2, SIFT, SNAP, and PANTHER.

Molecular Modeling and Homology
Modeling of the Mutated Genes By using MOE (Chemical Computing Group, MOE, v2014.0901, www.chemcomp.com) homology modeling was attempted in all seven genes that carried the missense variants identified in the present study. The selection of template crystal structures for homology modeling was based on the primary sequence identity and the crystal resolution. The MOE homology model method is separated into four main steps. First, a primary fragment geometry specification. Second the insertion and deletions task. The third step is the loop selection and the side-chain packing and the last step is the final model selection and refinement. The template selection was as follows: for PROP1 the 3A01 PDB file was used, for SRA1 the 4NBO PDB file was used, for RNF216 the 5l1V PDB file was used, for FGFR1 the 1CVS PDB file was used and for POLR3A the 5FJ8 PDB file was used. All models were handled, verified and visualized using the Drugster suite (23).

Model Optimization
Energy minimization for all four models was done in MOE initially using the Amber99 (24) force-field implemented into the same package, up to a root mean square deviation (RMSd) gradient of 0.0001 to remove the geometrical strain. The models were subsequently solvated with simple point charge (SPC) water using the truncated octahedron box extending to 7 Å from the model, and molecular dynamics was performed at 300 K, 1 atm with 2 fs step size for a total of 10 ns, using the NVT ensemble in a canonical environment (NVT stands for Number of atoms, Volume and Temperature that remain constant throughout the calculation). The results of the molecular dynamics simulation were collected into a database by MOE for further analysis.

Genetic Findings
All seven patients were sequenced by WES. The clinical, biochemical and genetic characteristics are summarized in Table 1. A total of nine variants were identified in genes that are known to be linked with the development of CHH/KS ( Table 2). Seven of these variants were novel and two were previously reported. The novel X-linked p.Gln82 * in the ANOS1 (KAL1) gene was found in patient 1, a 28-year-old CHH male with pubertal absence, cryptorchidism and micropenis (Table 1, Figure 1). The novel WDR11 p.Leu244Pro variant is probably pathogenic and is inherited in an autosomal dominant fashion (AD). This variant was identified in patient 2, a 72-year-old male with KS and associated clinical characteristics of anosmia, cryptorchidism and micropenis. Patient 2 first sought medical advice at the age of 40 years and since then remains a patient of our clinic. Molecular diagnosis was only possible 32 years later (Table 1, Figure 2). The previously reported AR, probably pathogenic p.Ile179Thr variant in the SRA1 gene was identified in heterozygosity in patient 3, a 19-year-old male with partial hypogonadism and upper limb defects (Table 1, Figure 3). In addition, the novel p.Asp792Asn in the RNF216 gene was also identified in heterozygosity in the same patient (Figures 3, 8,  Supplementary Figure 1). Evaluation by the ORVAL platform for digenicity predicted this novel variant to have a neutral effect (Supplementary Figure 6). However, familial segregation data and in silico structural models indicated a digenic mode of inheritance (Figures 3, 8, Supplementary Figure 5). Variants in the SRA1 and RNF216 genes have been associated with effects on the CHH phenotype. Thus, the presence in patient 3 of pathogenic, heterozygous variants in SRA1 and RNF216 genes could potentially be another example of digenic inheritance for the development of CHH (28). Patient 4, an 18-year-old male with CHH and the associated clinical characteristics of cryptorchidism and micropenis, was found to carry the novel AD p.Arg2400Trp variant in the CHD7 gene (Table 1, Figure 4). Patient 5, a 31-year-old male with CHH and associated clinical characteristics of cryptorchidism and micropenis, was found to carry the novel AD p.Pro186Ala in the FGFR1 gene (Table 1, Figure 5). The previously reported AD p.Arg822Cys also in the FGFR1 gene was found in patient 6, a 20year-old male with CHH and associated clinical characteristics   . The ORVAL platform was used to evaluate the pathogenicity of the variants regarding digenic inheritance, with respect to which they were shown to be positive (Supplementary Figure 6). Such variants when inherited in the AR form are usually associated with Pol III-related hypomyelinating leukodystrophies and not with CHH. Patients with Pol III-related leukodystrophies may have various clinical characteristics including ataxia, delayed dentition, hypomyelination, hypodontia, and hypogonadotropic hypogonadism. Patient 5, in addition to CHH, also developed a mild hypomyelinating leukodystrophy phenotype, which is likely associated with the heterozygous condition found in the POLR3A gene. It should also be noted that 3 of our patients were found to have other variants as well. The known p.Val103Ile variation of the MC4R gene, which has been linked to obesity, was found in heterozygocity in two of the male patients (the 18 and the 19-year-old males, Patients two and four). Both the patients were obese, with BMI above +2SDS and both developed insulin resistance. In addition the novel variant p.Arg112gGln in the PROP1 gene was found in the heterozygous state also in patient 4, an 18-year-old male (Figure 7, Supplementary Figure 4). Such variants have been reported to be associated with combined pituitary hormone deficiency. Our patient also had central hypothyroidism and he is currently on treatment with thyroxin.
The novel non-sense pathogenic variant p.Gln82 * in the ANOS1 gene, which was identified in a 28-year-old male, encodes a premature termination codon. It is expected to yield a truncated ANOS1 protein, missing the whey acidic protein (WAP)-like protease inhibitor domain and the four fibronectin type III (FN[III]) domains (Figure 1). The missense variants, identified by WES and confirmed by Sanger sequencing, were predicted to be deleterious by at least two prediction tools using the PredictSNP consensus classifier (22) ( Table 3). The PredictSNP consensus classifier evaluates the pathogenicity of a variant by using seven different prediction tools: MAPP, PhD-SNP, PolyPhen-1, PolyPhen-2, SIFT, SNAP, and PANTHER. Furthermore, the identified variants were absent from the population-specific data in an in-house database composed of 51 random samples of Cypriot origin. Familial segregation was available for two patients, patient 1 and patient 3 (Supplementary Figure 5, Table 2). For patient one the identified pathogenic variant followed the X-liked mode of inheritance and for patient 3, variants followed the digenic mode of inheritance (Supplementary Figure 5).

Conserved Protein Sequences Among Species
Protein alignment analyses of all identified pathogenic and probably pathogenic variants including the p.Gln82 * of the ANOS1 gene, the p.Leu244Pro of the WDR11 gene, the p.Asp792Asn of the RNF216 gene, the p.Arg2400Trp of the CHD7 gene, the p.Pro186Ala and p.Arg822Cys of the FGFR1 gene, the p.Ile179Thr of the SRA1 gene, the p.Arg561Gly of the POLR3A gene and the p.Arg112gGln of the PROP1 gene showed high amino acid conservation among different species (Figures 1B, 2B, 3B, 4B, 5B, 6B, 7B, 8C). Therefore, the newly and previously discovered variant sites of the above genes were probably located in a vital region of the coding genes that might affect the corresponding proteins mechanistically and/or structurally.
Effect of the Mutated PROP1, SRA1, FGFR1, POLR3A, and RNF216 at the Related Protein Structure The selection of template crystal structures for homology modeling was based on the primary sequence identity and crystal resolution models were obtained for the PROP1, SRA1, FGFR1, POLR3A, and RNF216. Unfortunately, no crystal resolution models were obtained for the CHD7 and WDR11 proteins.
The models for the SRA1, FGFR1, POLR3A, PROP1, and RNF216 were designed using means of homology modeling (Figures 3C, 5C, 6C, 7C, 8D). Each one of the recorded variants were induced in silico and each model was energetically optimized via Energy minimizations and Molecular Dynamics The wild type and variant models were subsequently superposed and the variant residues were inevitably superposed too The electrostatic surface was calculated for each wild type and variant set for the SRA1, POLR3A, and PROP1 models (Figures 3E, 6E, 7E) and finally the 2D interaction diagram was drawn for each pair to visualize all bonding and conformation changes induced upon by the variants (Figures 3F, 5E, 6F, 7F). Since the FGFR1 variant involved the replacement of a proline residue the conformational change in the 3D conformation of the beta sheet formation is showed in Figure 5F.
More specifically, the p.Arg112Gln variant in the PROP1 results in a significant physicochemical change in the 112 position. The bulkier and positively charged arginine residue is now replaced by a smaller polar and uncharged glutamine residue. Inevitably the interaction to the DNA molecule that was also included in the model is now lost, as well as a strong stabilizing H-bond to a nearby negatively charged aspartic acid residue. It is evident from the 2D interaction map that the mutant PROP1 has lost its potential to interact with DNA and that would unavoidably alter its molecular and cellular function (Figure 7). Moving on, the p.Ile179Thr variant of the SRA1 also changes significantly the physicochemical profile of the amino acid at 179 position. Isoleucine is a non-polar aliphatic residue, whereas threonine is a polar residue. This electrostatic change results in the establishment of two new bonds with two nearby leucine residues. This fixates the local three helical bundle conformation more strongly, thus making the 3D structural arrangement adjacent to the 179 position more compact and conformationally rigid compared to the wild type (Figure 3). Likewise, the p.Pro186Ala variant in the FGFR1 has a significant conformational impact. The proline residue is a very special amino acid that bears a cyclic side chain. This gives this residue a unique property of inducing a kink in the 3D conformational arrangement of the protein. Removing it, and replacing it with an alanine residue resulted in a 2.5 Angstrom shift outwards of the beta sheet it is located in and that consequently led to a 4.5 Angstrom shift of the neighboring antiparallel beta sheet formation too. Such significant changes in structural arrangements are bound to change the mechanics and function of the FGFR1 protein ( Figure 5). Lastly, the p.Arg561Gly in the POLR3A is a fine example of a variant change from a very bulky and positively charged residue to a small amino acid. From the 2D interaction diagram for this residue we can deduce that only the arginine amino acid is long enough to reach the nearby positively charged aspartic acid residue and to establish with the latter strong H-bonds. The glycine residue is not just small but it also misses the essential amino-groups that are required to establish H-bonding to the glutamic acid (Figure 6). The p.Asp792Asn variant according to the homology model of RNF216 is exposed to the solvent. It is located on a hairpin like loop, linking a beta-sheet and an a-helix in an antiparallel fashion ( Figure 8D). Therefore, there is high probability that this could be an interacting part of the RNF216 protein, judging from the rotamer of this residue, which is pointing outwards. In this direction, we modeled the electrostatic potential of the 792 residue position (Figure 8E), of the adjacent residues within a 7 Angstrom radius ( Figure 8G) and of the whole RNF216 protein (Figure 8F). It was found that the variation of aspartic acid to asparagine changes significantly the electrostatic potential and nature of the 792 position due to the extra NH2 moiety on the asparagine amino acid. That, coupled with the fact that all amino acids around the 792 position are negatively charged (red color- Figure 8G, left) is very significant as with the introduction of the asparagine residue a positively charged group is now introduced. Taking the abovementioned facts into consideration, we propose that the p.Asp792Asn variant is very significant as the positively charged moiety that is introduced could disrupt the conserved negatively charged region of RNF216, thus leading to a considerable change in the physicochemical and electrostatic profile of that domain that would inevitably affect its binding/interaction potential and would probably change its functional properties (e.g., loss of recognition and even interaction).

DISCUSSION
The present study investigated by high-throughput whole exome sequencing the genetic impact in patients with CHH. The seven patients of Cypriot origin with CHH/KS were identified with variants in genes linked with this phenotype: ANOS1, SRA1,  CHD7, WDR11, FGFR1, RNF216, and POLR3A. A total of seven novel and two rare previously reported variants were identified in the patients of the current study and were found as novel or very rare in the ExAC population database (29). All these variants were also predicted to be pathogenic by at least two computational programs (22,(30)(31)(32)(33)(34)(35)(36). Our results once more confirmed the genetic complexity of CHH and the roles that exemplify a series of pleiotropic genes during development (13,37,38).  More specifically, we identified the novel X-linked hemizygous truncated p.Gln82 * pathogenic variant in the ANOS1 gene in a 28-year-old male with CHH ( Table 1). Patients with sporadic KS/CHH due to ANOS1 gene defects have been correlated with the phenotype of right renal agenesis/dysgenesis, thus provide evidence for the X-linked mode of inheritance and offering the opportunity for genetic counseling (12,(39)(40)(41)(42)(43)(44). Approximately 10-20% of males with KS carry ANOS1 variants or intragenic microdeletions (38). The majority of the X-linked KS variants cause alteration of splicing, frameshift or stop codons leading to synthesis of truncated anosmin (45). Nonsense variants and full deletions in the ANOS1 gene are the most pathogenic and lead to a truncated and absent anosmin protein, respectively (37). Two of our normosmic male adult patients have been identified with the novel AD p.Pro186Ala and the previously reported p.Arg822Cys (46) variants in the FGFR1 gene. Both of these patients were also characterized by delayed puberty during adolescence and later CHH. In a similar fashion with ANOS1, the expression of the FGFR1 gene is also generated in the apparent olfactory bulbs and loss-of-function variants cause a form of KS with autosomal dominant inheritance (4-24, 28-50). FGFR1 is a cell surface membrane receptor that possesses tyrosine kinase activity and mediates fibroblast growth factor signaling (51). Patients with variants in the FGFR1 gene present also various congenital anomalies that are not associated with the reproductive system and are often associated with kidney and tooth differentiation, ear and palate morphogenesis and the development of cortico-spinal axonal tracts (52). Notably, patient five, the 31-year-old male patient with CHH identified with the novel FGFR1 p.Pro186Ala also shared in heterozygosity the novel POLR3A p.Arg561Gly missense variant. This finding adds to the already known spectrum of phenotypes resulting from POLR3A and POLR3B variants. POLR3A and POLR3B can be also associated to neurological or dental anomalies and isolated hypogonadotropic hypogonadism (53). Patient 4 in addition to the novel AD CHD7 p.Arg2400Trp variant also carried the known MC4R p.Val103Ile variant implicated in BMI and the novel p.Arg112gGln in the PROP1 gene. Various studies have described PROP1 gene variants as responsible for causing combined pituitary hormone deficiency (54)(55)(56). Heterozygous autosomal dominant loss-of-function variants in the CHD7 gene are the major causal factor of CHARGE syndrome (57,58), in addition to the fact that CHD7 variants have been also been reported in patients with isolated CHH (59)(60)(61). Several reports have also linked PROP1 variants with gonadotroph function that progressively declines and clinically patients with such variants may manifest a shortage of pubertal development, i.e., failure to enter or complete puberty (62,63). There are several reports of spontaneous puberty with a posterior decline of gonadotrophic function that have been linked to p.Arg112Ter, p.Arg120Cys, p.Phe88Ser, and c.150delA PROP1 variants (64)(65)(66)(67). Since the PROP1 gene is involved in the anterior pituitary, cell lineage specification variants could behave as an additive factor in the development of CHH when co-inherited with variants from genes involved in normal gonadotroph function. Such could be the case with the 18-yearold CHH patient of the present study identified with the novel CHD7 p.Arg2400Trp and the novel p.Arg112gGln variant in the PROP1 gene.
In the present study, a 72-year-old anosmic KS patient originally sought medical advice at the age of 40 in our clinic. Since then, he remains a patient of our clinic and at the age of 72-years he was identified with the novel AD p.Leu244Pro in the WDR11 gene. WDR11 has been implicated in CHH and KS, human developmental genetic disorders defined by delayed puberty and infertility (68,69). Several reports in CHH patients with and without anosmia identified in heterozygosity variants in the WDR11 gene (68). WDR11 is expressed in several adult organs including the brain and the gonads. Comprehensive analysis of the mouse brain displayed WDR11 expression in the GnRH neuronal migratory location including nasal cavity and cribriform plate area in E12.5 mouse embryo as well as the median eminence in the adult brain, showing co-localization with GnRH. Furthermore, WDR11 is expressed all over the developing and adult olfactory bulb (OB) and its WD domains are important for β-propeller formation and protein-protein interaction (70). In addition, WDR11 interacts with EMX1, a homeodomain transcription factor involved in the development of olfactory neurons, and missense variants diminish or eliminate this interaction (68). Therefore, it is highly likely that the impaired pubertal development in these patients results from a deficiency of productive WDR11 protein interaction.
Interestingly, two out of the seven CHH patients in our cohort, a 30-year-old female (Patient 7) and a 19-year-old (Patient 3) male were both identified with variants in the SRA1 gene. More specifically, the 30-year-old female carried in homozygosity the previously reported p.Ile179Thr variant in the SRA1 gene (25). The 19-year-old male also carried this same SRA1 variant in heterozygosity together with the novel p.Asp792Asn variant in the RNF216 gene. As reported by Kotan et al. (25), the variant p.Ile179Thr was reported only once in one independent Turkish family with IHH/delayed puberty and its severity was supported by functional studies. Using a mutant SRA1 construct, reduced co-activation of ligand-dependent activity of the estrogen receptor alpha was demonstrated (25).
The variant p.Ile179Thr was not found in 51 Cypriots used as controls for the purposes of the present study and was reported with an allele frequency of 0.00081 in GnomAD v2.1.1. Therefore, most likely, a hot spot exists for this specific variant in the greater Eastern Mediterranean region, suggesting a founder effect phenomenon, which has been also seen for other rare endocrine conditions in this area (71).
The SRA1 is a steroid receptor RNA activator that has been shown to positively regulate the activity of the androgen receptor and the estrogen receptor (72,73). In recent years only a few studies linked the SRA1 gene as responsible for causing CHH when patients inherit pathogenic variants in the AR form (25,26).
The concepts of incomplete penetrance and variable expressivity have been notified in such cases as the 19-year-old patient where digenic variants are observed. Digenic variants account for variable phenotypes in idiopathic hypogonadotropic hypogonadism and other disorders and several recent and older reports identified such conditions (17,40,74,75). The existence of digenic and oligogenic inheritance in CHH is quite common, with about 20% of CHH cases reported to share at least two causative variants that could result in a disease phenotype (1,15). The most appropriate way to examine the possibility of low penetrance and variable expressivity of CHH genes is by concurrently carrying out targeted genetic analyses, or preferably by performing WES on the probands and available relatives so as to establish digenic or oligogenic transmission.

REVIEW OF THE LITERATURE
Congenital hypogonadotropic hypogonadism (CHH) is a rare disorder of sexual maturation characterized by GnRH deficiency with low sex steroid levels associated with low levels of LH and FSH. CHH may be caused by variants in numerous genes and recent studies shed light on the complexity of CHH genetics (6,15,76). Over the recent few years, genetic evaluation of patients with inherited diseases, including CHH, has increasingly utilized massive parallel sequencing by next-generation sequencing (NGS), that allows the concurrent investigation of thousands of genes (1,15,77). At this scale of analysis, NGS is inexpensive and rapid compared to the traditional Sanger sequencing and is increasingly being used in medical practice. NGS has certainly facilitated CHH genetic diagnoses and aided healthcare professionals to provide reliable and informed genetic counseling for patients with CHH. The crucial challenge regarding NGS concerns the identification of true oligogenism in circumstances involving several rare variants which do not have a clear phenotypic effect and are identified by coincidence. Such a challenge also concerns the identification of genes underlying CHH pathogenesis and which are likewise reported to act in an oligogenic context (78). Since the discovery of ANOS1 (79), more than sixty genes have been reported to underlie CHH and were previously considered to be inherited in the AD form (6). Herein, we review six of these genes: ANOS1, FGFR1, CHD7, WDR11, RNF216, and POLR3A, since novel variants in these genes have been identified in our cohort of patients under investigation.

ANOS1 (KAL1) Gene Variants Causing X-Linked Recessive KS/CHH
ANOS1 was the first gene linked to Kallmann syndrome (KS) and since the early nineties when the first reports demonstrated variants with an X-linked mode of inheritance (80-83), many others followed throughout the years (6,(84)(85)(86)(87)(88). KS occurs more frequently in males than in females, with an estimated prevalence of 1 in 30,000 males and 1 in 120,000 females (12). Patients with KS associated with ANOS1 pathogenic variants usually exhibit anosmia accompanied with CHH (12,14,85,86). Fewer patients with pathogenic variants in ANOS1 are either anosmic or hyposmic and have been reported to exhibit other signs, such as mirror movements and renal agenesis, but they do not always co-segregate with the variant recognized in a given family (85,(89)(90)(91). According to the Human Gene Mutation Database (http://www.hgmd.cf.ac.uk/ac/index.php) more than 150 ANOS1 pathogenic variants have been reported as the causative factor in KS patients. Most of these pathogenic variants mainly consist of nucleotide deletions or insertions and to a lesser extend of variants that involve amino acid missense substitutions (88,92,93). The ANOS1 gene encodes anosmin, a protein which plays a significant role in the embryogenesis of brain, kidneys, respiratory and digestive systems (92,94). Anosmin, as an extracellular matrix protein binds to the cell membrane and stimulates the development of the olfactory system and behaves as an axonal guidance for the GnRH neurons, the olfactory cells and the Purkinje cerebellum neurons (95). Monogenic loss-offunction pathogenic variants in ANOS1 gene have been estimated to account for 4-10% of KS/CHH cases and has been principally studied in many reported cohorts (12, 41-44, 88, 96-99).
Regarding the reproductive phenotype, male KS patients with ANOS1 variants display a complete penetrance of CHH and their pre-and postnatal gonadotropin deficit is severe with a high frequency of micropenis, cryptorchidism and complete absence of gonadal development (15,16).

FGFR1 Gene Variants Causing KS/CHH
The presence of variants in the FGFR1 gene is another important cause of KS and was the first gene to be identified as an AD form of the disease (49,100). More than 140 loss-of-function mutations in the FGFR1 gene have been reported with missense, non-sense and frameshift defects being the most frequent (101). Less frequent autosomal gene deletions have also been reported in patients with CHH and KS (102,103). FGFR1 is considered to be a pleiotropic gene that can display different roles during development and variants found in it can cause CHH with or without anosmia (49,100). Genotype-phenotype correlations in patients with AD variants in the FGFR1 gene demonstrated some clinical features linked with KS, such as loss of nasal cartilage, hearing deficit and anomalies of the limbs (6,93,101). The function of FGFR1 in the normal development of the olfactory bulb proposes the link of anosmia with GnRH deficiency in the FGFR1-mutated patients (104). Phenotype analysis proposes that FGFR1 elaborates in the normal migration of GnRH fetal neurons, but this is not entirely clear-cut as a considerable proportion of FGFR1-mutated patients have normosmic GnRH deficiency (15). Regarding the reproductive phenotype of male patients with FGFR1 variants, the penetrance of CHH and GnRH deficiency is variable and ranges from profound to partial puberty and even to reversal (1,99,105,106).
Several groups have reported patients harboring FGFR1 variants linked to non-reproductive signs. Patients with FGFR1 mutations have been reported to suffer from health conditions such as 8p11 myeloproliferative syndrome (107,108), encephalocraniocutaneous lipomatosis (109,110), Hartsfield syndrome, a rare condition characterized by holoprosencephaly, which is an abnormality of brain development (111,112), osteoglophonic dysplasia, a condition characterized by abnormal bone growth that leads to craniofacial abnormalities and dwarfism (113,114) and Pfeiffer syndrome, which is characterized by craniosynostosis (115,116). Somatic pathogenic variants involving the FGFR1 gene have also been reported in several types of cancers, including the lung, breast, esophagous, oral cavity and brain tumors (101,(117)(118)(119). Taking into consideration the genotypic and phenotypic heterogeneity that is observed in patients with FGFR1 variants and the fact that their prevalence is not clearly established makes genetic counseling rather complicated.

CHD7 Gene Variants Causing CHARGE Syndrome and CHH/KS
CHD7 is the gene that codes for the chromodomain helicase DNA binding protein 7 and variants found in the AD form were first reported as the genetic cause in a series of patients with CHARGE (coloboma, heart defect, atresia choanae, growth retardation, genital abnormality, and ear abnormality) syndrome (57,120). CHARGE syndrome occurs in approximately 1 in 8,500 to 10,000 new-borns and up-to-date more than 600 CHD7 variants in the AD form have been associated with the disorder (57,(121)(122)(123)(124). Several other studies of families carrying CHD7 mutations in the AD form also demonstrated a broad phenotypic variability and linked more than 50 of them with KS and congenital hypogonadotropic hypogonadism (59,60,69,125). It has been estimated that inherited and de novo CHD7 mutations account for 5-10 percent of all cases of KS and an accountable number of these patients exhibit mild form features of CHARGE syndrome, such as abnormally shaped ears, hearing loss, hare lip/cleft palate and cardiac abnormalities (60,61,69,121,126).

WDR11 Gene Variants Causing CHH/KS
WDR11 is a member of the WD-repeat containing protein family and comprises of twelve conserved domains of approximately 40 amino acids (68). The WDR11 gene is located in the chromosome 10q25-26 region and is expressed in various human organs including the brain, ear, lung, heart, kidney and the gonads (70). WDR11 is a scaffolding protein that is involved in multiple of cellular proceedings, including cell cycle progression, signal transduction, apoptosis and gene regulation (70). Kim et al. first reported that when mutated WDR11 is linked with idiopathic HH and KS (68). Since the initial report by Kim et al. (68), a few others followed and linked the WDR11 gene with different pathogenic variants in male patients without anosmia and CHH (127,128). Recently, the WDR11 gene has also been shown to be involved in the Hedgehog (Hh) signaling pathway which is important for the normal ciliogenesis and when mutated can be the causal factor of KS and HH (68,70). Another recent report by Sutani et al. (129) linked WDR11 as another causative gene for coloboma, cardiac anomaly and growth retardation in the 10q26 deletion syndrome.

RNF216 Gene Variants Linked to Gordon Holmes Syndrome
The RNF216 protein is a cytoplasmic protein which interacts with the serine/threonine protein kinase i.e., the receptor-interacting protein (RIP). Particular zinc finger domains of the RNF216 protein are necessary for its interaction with RIP and for the inhibition of TNF-and IL1-induced NF-kappa B activation pathways (130,131). Additionally, the RNF216 protein plays a role in the ubiquitin-proteasome system for the break-down and degradation of unwanted proteins. Specifically, this protein functions as an E3 ubiquitin ligase (132). Variants in the RNF216 gene have been linked with hypogonadotropic hypogonadism, ataxia and dementia (28). More explicitly digenic homozygous variants in RNF216 and OTUD4, which encode a ubiquitin E3 ligase and a deubiquitinase, respectively, were identified in three affected siblings in a consanguineous family (28). Several other recent studies also reported variants in the RNF216 gene as a result of consanguinity to cause Gordon Holmes syndrome, a rare disorder characterized by diminished production of hormones leading to hypogonadotropic hypogonadism and difficulty in the coordination of movements i.e., cerebellar ataxia (133)(134)(135)(136). These recent findings regarding the RNF216 gene associate the disorderly ubiquitination to neurodegeneration and reproductive dysfunction in combination with functional studies to reveal specific genetic interactions that cause disease.

POLR3A Gene Variants Associated With Hypomyelinating Leukodystrophy and HH
The POLR3A gene provides instructions for the production of the largest subunit of RNA polymerase III which is the enzyme involved in the RNA synthesis (137). The gene is located in chromosome 10q22.3 and variants inherited in the AR form have been initially reported in French-Canadian families with hypomyelinating leukodystrophy (138). Interestingly, these families were mapped to the same locus as leukodystrophy with oligodontia and demonstrated clinical and radiological overlap with patients with hypomyelination, hypodontia and hypogonadotropic hypogonadism syndrome (138). Several other recent studies that followed also reported variants in the POLR3A gene as being responsible for causing hypomyelination, hypodontia and hypogonadotropic hypogonadism, thus establishing a series of POLR3A gene variants to be associated with polymerase III-related leukodystrophy (139)(140)(141)(142)(143)(144). It is estimated that 30-40% of patients with leukodystrophies remain without a molecular diagnosis (138,141). The existence of mild and overlapping hypomyelinating leukodystrophy phenotypes could be attributed to heterozygous variants found in the POLR3A gene as a result of an abnormal enzymatic function of the RNA Pol III catalytic subunit. The role of heterozygosity in POLR3A in the overall pathogenesis of CHH is not well-established, and the possibility of a synergistic effect between these variants and variants identified in other genes cannot be excluded. Additionally, POLR3A gene could also be speculated to be a phenocopy gene due to the observed variability of phenotypes, therefore, patients and family members identified with mutations in this gene should be re-evaluated for understated and previously unrecognized clinical signs.

CONCLUSION
GnRH deficiency has been recognized both clinically and genetically as a heterogeneous disease with a range of different reproductive phenotypes including of congenital GnRH deficiency with anosmia (KS) and congenital GnRH deficiency with normal olfaction (normosmic CHH). The present study/review discusses the involvement of known and novel variants in patients with CHH/KS and adds up to the ontogeny of GnRH deficiency.
Moreover, this study provides new genetic findings and reinforces the significance of the use of NGS technology for the accurate molecular diagnosis and treatment of this rare condition.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by The Cyprus National Ethics Committee. Written informed consent to participate in this study was provided by the participants' legal guardian/next of kin. Written, informed consent was obtained from all seven adult individuals that participate in the study for the publication of any potentially identifiable images or data included in this article.

AUTHOR CONTRIBUTIONS
VN, PF, MT, NS, and LP contributed to the conception, design, and interpretation of the project. MS, FC, AO, GS, and DV contributed to the experimental part and data analysis of the project. NN, GT, and CM contributed to drafting or revising intellectual content of the manuscript. VN, NS, and LP had primary responsibility for final content. All authors read and approved the final manuscript. Black arrows indicate the disease-causing position of the digenic combination on the plot. Y-axis, Support Score: percentage of individual predictors agreeing on the disease-causing class of the digenic combination. X-axis, Classification Score: median probability among all predictors that the digenic combination is disease-causing. The higher the scores, the more confident the predictor is for the disease-causing class.
Supplementary Table 1 | Genes with biological involvement in the GnRH neuronal system and CHH selected for variant analysis of the WES. For each gene the corresponding genetic and phenotypic information is indicated. The mean target coverage of the target region for each gene is designated in a separate column.