EAP1 regulation of GnRH promoter activity is important for human pubertal timing

Abstract The initiation of puberty is orchestrated by an augmentation of gonadotropin-releasing hormone (GnRH) secretion from a few thousand hypothalamic neurons. Recent findings have indicated that the neuroendocrine control of puberty may be regulated by a hierarchically organized network of transcriptional factors acting upstream of GnRH. These include enhanced at puberty 1 (EAP1), which contributes to the initiation of female puberty through transactivation of the GnRH promoter. However, no EAP1 mutations have been found in humans with disorders of pubertal timing. We performed whole-exome sequencing in 67 probands and 93 relatives from a large cohort of familial self-limited delayed puberty (DP). Variants were analyzed for rare, potentially pathogenic variants enriched in case versus controls and relevant to the biological control of puberty. We identified one in-frame deletion (Ala221del) and one rare missense variant (Asn770His) in EAP1 in two unrelated families; these variants were highly conserved and potentially pathogenic. Expression studies revealed Eap1 mRNA abundance in peri-pubertal mouse hypothalamus. EAP1 binding to the GnRH1 promoter increased in monkey hypothalamus at the onset of puberty as determined by chromatin immunoprecipitation. Using a luciferase reporter assay, EAP1 mutants showed a reduced ability to trans-activate the GnRH promoter compared to wild-type EAP1, due to reduced protein levels caused by the Ala221del mutation and subcellular mislocation caused by the Asn770His mutation, as revealed by western blot and immunofluorescence, respectively. In conclusion, we have identified the first EAP1 mutations leading to reduced GnRH transcriptional activity resulting in a phenotype of self-limited DP.


Introduction
Puberty represents the remarkable transition from childhood to adult life with the attainment of reproduction and adult stature. The onset of puberty is triggered by the reactivation of the hypothalamic-pituitary-gonadal (HPG) axis through augmentation of the pulsatile release of gonadotropin releasing hormone (GnRH) from the hypothalamus, which in turn causes an increased pulsatile gonadotropin release from the anterior pituitary. Pubertal onset shows a large variability among populations, due to ethnicity, environment, nutrition and stress, (1)(2)(3) although previous studies indicate that 60-80% of the variance in pubertal timing can be explained by heritable factors (4,5). However, the precise mechanisms underlying the regulation of puberty onset are still unknown.
In 2000, Rampazzo et al. (6) identified an intronless gene containing a zinc (Zn) finger domain, mapped to 14q24.3 and provisionally named C14ORF4. Following functional annotation of C14ORF4 that indicated its pivotal role in pubertal timing, it was renamed enhanced at puberty 1 (EAP1) (7) and later was assigned the name of interferon regulatory factor 2 binding proteinlike (IRF2BPL). EAP1 mRNA has been shown to increase in the hypothalamus of rats and non-human primates at the time of puberty, and EAP1 deficiency led to delayed puberty (DP) and disrupted estrous cyclicity in both rodents (7) and nonhuman primates (8). Eap1 codes for a nuclear transcription factor, characterized by a dual transcriptional activity; it both transactivates the GnRH promoter, which facilitates GnRH secretion, and inhibits the preproenkephalin promoter, which represses GnRH secretion. Therefore, Eap1 transcriptional activity facilitates the initiation of female puberty, in a manner that is independent of hypothalamic Kiss1 expression (9). Notably, a study performed in rats showed that hypothalamic expression of Eap1 is not directly regulated by ovarian steroids, as its expression in the peri-pubertal female hypothalamus changes even in the absence of ovaries (10). Despite this seemingly important role, no EAP1 mutations have yet been identified in humans with pubertal disorders. DP is defined as the absence of testicular enlargement in boys or breast development in girls at an age that is 2 to 2.5 standard deviations (SDs) later than the population mean (11). In absence of any identifiable cause, DP usually resolves by the age of 18 years, and in this case is referred to as self-limited, or constitutional, DP (12). Self-limited DP is commonly familial and segregates with an autosomal dominant inheritance patterns in >70% of families, indicating a strong genetic basis of the trait. Such an inheritance pattern suggests that DP has a monogenic or oligogenic background, although very few underlying genes have been discovered (13). Our large, well-phenotyped Finnish DP cohort (14) offers a prodigious tool to investigate potential pathogenic variants causing self-limited DP, likely enriched in our cohort. Indeed, comparable next generation sequencing strategies in patients from this cohort have unveiled roles for novel genes in the pathogenesis of self-limited DP (15,16).
In this study we have identified, for the first time, two EAP1 mutations leading to reduced GnRH transcriptional activity resulting in the phenotype of self-limited DP.

Exome sequencing of families with self-limited DP identifies variants in EAP1
Whole-exome sequencing (WES) of 67 informative families from our large cohort with self-limited DP returned 6 952 773 variants (Fig. 1A). These variants were filtered through our in-house pipeline to identify rare, predicted deleterious mutations, segregating in multiple families and with potential biological relevance. The 28 top candidate genes identified then underwent targeted resequencing in additional 42 families from the same cohort (178 individuals with DP and 110 controls), and the filtered results were analyzed by applying statistical thresholds for enrichment of rare, pathogenic variants in our cohort via whole-gene rare variant burden testing (RVBT) with multiple comparison adjustment (17). Four genes passed the P < 0.025 threshold on RVBT, and potentially pathogenic variants in these genes were further analyzed to determine their presence in controls from our cohort and for segregation with trait (Fig. 1A).
The candidate gene, EAP1 (ENSG00000119669; gene identification number 64207; synonyms: IRF2BPL and Chromosome 14 Open Reading Frame 4 [c]), was identified after RVBT (adjusted P = 1.99E-05). EAP1 has been highlighted as an important potential candidate for the hierarchical regulation of pubertal timing through systems biology approaches and animal models (18). EAP1 mRNA levels and protein expression are seen to increase in the hypothalamus of primates and rodents at the time of pubertal onset (7).
Two identified variants in EAP1 are rare, highly conserved and potentially damaging to protein function Eight rare and potentially pathogenic variants were identified in EAP1 from whole and targeted exome sequencing. Five of these variants were discarded in our post-sequencing analysis, as they were present in multiple controls from our cohort. Three variants of interest in EAP1 were initially identified in four pedigrees from the cohort, but only two were found to segregate within families after Sanger sequencing of all family members. These two variants [NM 024496.3: c.2308A>C (rs760847179) p.Asn770His Chr: 14:77491828 and NM 024496.3: c.661 663delGCG p.Ala221del Chr: 14:77493473] were each identified in one proband from the cohort and their affected family members. The two probands from these pedigrees did not carry any other predicted pathogenic variants in known GnRH deficiency or gonadotropin deficiency-causing genes (16). The two EAP1 variants, p.Asn770His and p.Ala221del, are both rare with a minor allele frequency (MAF) of <0.5% in population databases. Both variants affect amino acids that are highly conserved among homologues, as revealed by PhyloP score, and multiple sequence alignment (Fig. 1B). The p.Asn770His missense variant is predicted to be deleterious to protein function by prediction tools and CADD score, as the affected amino acid residue resides in a long loop region, which is C-terminal to the C3HC4 ring domain (residues 715-762). This region shows evolutionary conservation among different species (Fig. 1C). Although the amino acid substitution (p.Asn770His) does not affect the residues directly involved in Zn binding and it is not predicted to reduce protein stability, it introduces a positively charged residue in place of a neutral one at a potential proteinprotein interaction site and is thus likely to be disruptive. The p.Ala221del is a novel deletion located within a predicted disordered and glycine-proline-rich region of the protein, but unfortunately no structural information could be obtained.

Family members with classical self-limited DP from two pedigrees carry heterozygous EAP1 variants inherited in an autosomal dominant pattern
Segregation with a clear autosomal dominant pattern of inheritance was seen in both the identified families, with the Filtering strategy to identify EAP1 as a candidate gene for self-limited DP, prediction of pathogenicity and conservation across species for EAP1 variants identified, and structural model of the Zn finger domain of EAP1. (A) WES was performed on DNA extracted from peripheral blood leukocytes of 160 individuals from our cohort (67 DP probands, 58 DP relatives and 35 controls). The exome sequences were aligned to the UCSC hg19 reference genome. Picard tools and the genome analysis toolkit were used to mark PCR duplicates, realign around indels, recalibrate quality scores and call variants. Variants were filtered for potential causal variants using filters for quality control, predicted functional annotation, MAF, case-control analysis, variants in multiple families and biological relevance. Targeted exome sequencing using a Fluidigm array of 28 candidate genes identified post-filtering was then performed in a further 42 families from the same cohort (288 individuals, 178 with DP and 110 controls). Variants post-targeted resequencing were filtered using the same criteria as the WES data. RVBT was performed with a multiple comparison adjustment applied post-hoc (17). Screening of 100 further cohort controls was via conventional Sanger sequencing. * Data published (15,16). * * Excluded due to the presence of variants in multiple controls. (B) Minor allele frequencies for ExAC Finnish population (accessed February 2018), conservation and pathogenicity scores (SIFT (37); Polyphen2 (36); MutationTaster (46)). Multiple sequence alignment (msa) was generated using MutationTaster (46). The p.Ala221 residue is highly conserved among different species, PhyloP score 1.801. The p.Asn770His is highly conserved among different species and the PhyloP score is 4.523. (C) The mutant 3D structure of the Zn finger domain of EAP1 is presented as a cartoon. The Zn atoms are presented as magenta spheres, the conserved C3HC4 residues (which bind Zn atoms) are presented in green and the mutant histidine (H) at position 770 is presented in red. The mutant His770, shown in red, is located on the surface of EAP1 and may be part of a protein-protein interaction site. The position of invariable residues C3HC4 is indicated with green arrows, whereas the position of N770 is indicated with a red arrow. The amino acid numbering and secondary structure is presented above the msa. C indicates cysteine. affected individuals carrying heterozygous changes in EAP1 variants ( Fig. 2A and B). The affected individuals from these two families have classical clinical and biochemical features of self-limited DP, with delayed onset of Tanner stage II and delayed peak height velocity (PHV) ( Fig. 2A , A , B and B ). The proband from family A presented at the age of 15.7 years with delayed pubertal development. His father had had similar delay in puberty onset. At initial evaluation the proband had testes volumes of 4 ml and circulating testosterone concentration at early pubertal level (5.7 nmol/l). The proband from family B presented at the age of 15.3 years with prepubertal testes volumes of 3 ml bilaterally. Both probands had markedly delayed bone age at presentation and during follow-up they had spontaneous pubertal development without testosterone therapy excluding idiopathic hypogonadotropic hypogonadism.

Eap1 is expressed within key regions within the mouse adult hypothalamus
We performed in situ hybridization on peri-pubertal hypothalamus of male and female mice. An abundance of Eap1 mRNA was detected in the ventromedial (VMH), paraventricular (Pa) and arcuate (Arc) nuclei of male (Fig. 3A) and female (Fig. 3B) mice. Eap1 expression was also investigated concurrently with the detection of GnRH neurons using immunohistochemistry. GnRH neuronal axons were mainly detected at the level of the median eminence (ME; Fig. 3A), whereas GnRH neuron bodies were interspersed in the median preoptic area in a positive Eap1 milieu ( Fig. 3C; Fig. 3D image shows GnRH neuron bodies at high magnification). Eap1 expression specificity was tested with a sense probe (Fig. 3E), which resulted in no detectable staining. Hence, Eap1 displayed strong expression in specific  (15). TH means target height based on mid-parental heights. Green dots connected by continuous black lines indicate bone age at the corresponding chronological age (blue dots), as estimated by the Greulich and Pyle method. hypothalamic subregions of both female and male peri-pubertal mice. These results in mice are in keeping with Eap1 expression in rat and non-human primate hypothalami (7). Previous evidence showed GnRH neurons expressing Eap1 (7).

EAP1 binding to the GnRH promoter increases at puberty
Results of a chromatin immunoprecipitation (ChIP) experiment indicated that EAP1 binds to GnRH1 promoter (Fig. 4) and that association of EAP1 to the GnRH1 promoter region in the medial basal hypothalamus (MBH) of female monkeys is increased at the late juvenile (LJ), as compared to the early juvenile (EJ) period, i.e. at the initiation of puberty. No change in EAP1 binding was detected in intron 2 of the GnRH gene. ChIP also performed using an antibody to beta-galactosidase (not expressed in brain tissue) did not show any appreciable signal in either EJ or LJ animals.

EAP1 variants significantly impair GnRH promoter activity in a dose-dependent manner
Previous results (7), together with evidence from our ChIP analysis, demonstrate that EAP1 transcriptional activity is promoter specific. We therefore employed a promoter assay in order to study GnRH promoter activity when trans-activated by EAP1 WT or mutated constructs (Fig. 5A). In a HEK293T cell line, the trans-activating strength of EAP1 on the human GnRH promoter was significantly reduced by the p.Ala221del in frame deletion variant (adjusted P = 0.0374) and highly significantly reduced by the p.Asn770His missense variant (adjusted P = 0.0003), as compared to the WT protein. The deletion of the RING finger domain (RINGdel), which is required for EAP1 transcriptional activity (7), impairs the ability of EAP1 to trans-activate the GnRH promoter (adjusted P = 0.00021) and was used as a control (Fig. 5A). In addition, as these are heterozygous mutations in our human patients, we co-transfected in equal amount the EAP1 WT expression vector with either the EAP1 Ala221del or with EAP1 Asn770His vectors, and we observed a dose-dependent reduction of EAP1 transcriptional activity (WT versus WT+Ala221del, adjusted P=0.0498; WT versus WT+Asn770His, adjusted P=0.0305).

EAP1 mutant proteins levels are altered compared to wild-type protein
To further understand the mechanism by which these mutations affect protein function, we examined protein expression in a HEK293T cell line, which has no endogenous EAP1 expression (Fig. 5B). EAP1 WT and mutated proteins were expressed (Fig. 5B) at the expected molecular weight of 90 kDa, while the RINGdel protein was expressed at 80 kDa, as expected due to the deletion of amino acids 715-762 (7). Densitometry analysis of the immunoblot (Fig. 5C), revealed that the Ala221del mutant protein levels, normalized to the housekeeping gene glyceraldehyde-3-phosphate dehydrogenase (GAPDH), are significantly reduced compared to WT (adjusted P = 0.0062), thus explaining the reduced activity of this mutant protein. However, the Asn770His mutant protein levels are increased significantly (adjusted P = 0.0322).

The EAP1 Asn770His mutant protein is mislocated in the cytoplasm
EAP1 is a transcriptional factor that localizes to the nucleus. To determine whether the Ala221del or Asn770His mutations affect the ability of EAP1 to reach the cell nucleus, we performed an immunocytoflourescence analysis (Fig. 6) on HEK293T cells transiently transfected with EAP1 WT, Ala221del or Asn770His mutant vectors. The results showed that both the EAP1 WT protein and the Ala221del mutant are expressed in the nucleus, as evidenced by co-expression with the nuclear marker 4 ,6diamidino-2-phenylindole (DAPI; Fig. 6C and F). In contrast, the Asn770His mutant protein failed to reach the nucleus, remaining in the cytoplasm (Fig. 6I), hence providing a compelling explanation for the reduced functional activity of this mutant.

Discussion
The central control of pubertal onset, after a period of juvenile quiescence, is mediated by a resurgence of the GnRH pulse generator, with a profound increase in the activity of the HPG axis. Although kisspeptin secretion from KNDy neurons in the Arc nucleus is one of the most important stimulatory elements of this neural network, kisspeptin has not been demonstrated . EAP1 binding to the GnRH promoter increases at puberty of female rhesus monkeys. An increased association of EAP1 to the rhesus GnRH1 promoter region in the MBH of female monkeys between the EJ and LJ periods is demonstrated. No change in EAP1 binding is detected in intron 2 of the GnRH gene. Inp means Input DNA and bGal means ChIP performed using an antibody to beta-galactosidase (a protein not present in the brain), serves as a negative control. Experiment was performed once, due to the limitation of availability of primate tissue.
as the ultimate controller of the release of the puberty brake but instead is likely to act as a conduit for upstream regulators (19). A growing body of evidence has demonstrated that no single gene alone is likely to be responsible for hypothalamic control of puberty. Instead, a sophisticated network of transcriptional factors, hierarchically organized, has been proposed as the machinery governing the balance between inhibitory and excitatory upstream inputs on the GnRH system (20). Transcriptional repression is a fundamental component of the regulation of mammalian gene expression, and transcriptional repressors containing Zn finger motifs, which recognize specific DNA sequences in regulatory regions of the genome, are particularly appealing candidates to have major roles in this governing network (21).
In this study we have identified through whole and targeted exome sequencing two deleterious mutations in EAP1 in pedigrees with self-limited DP. Eap1 had been proposed in several animal studies to have an important role in regulating the time of puberty onset; however, mutations in this gene have not, to our knowledge, been previously identified in conditions of abnormal pubertal timing in humans. EAP1 mutations were significantly enriched in our cohort of patients with self-limited DP as compared to control populations, and these mutations were inherited in the expected autosomal dominant pattern seen in this condition (14,22). The affected members of these families displayed typical self-limited DP with late onset of puberty but full adult development by the age of 18 years. None of the affected individuals had neurological or other associated phenotypic abnormalities (23). Complementing previous studies in rats and non-human primates (7), we have shown that Eap1 is also abundantly expressed in the peri-pubertal mouse hypothalamus. Previously published work provides evidence that GnRH neurons express Eap1 (7). Moreover, we have demonstrated that EAP1 binding to the GnRH1 promoter increases at the onset of puberty in female monkeys. In a previous rodent model siRNAmediated Eap1 knockdown was shown to cause delayed vaginal opening in female rats in addition to disrupted estrous cyclicity; reduced plasma luteinising hormone (LH), follicular stimulating hormone (FSH) and estradiol levels; and delayed growth of ovarian follicles (7). In non-human primates knockdown of EAP1 expression in the Arc lead to cessation of menstrual cyclicity, and a single-nucleotide polymorphism in the 5 -flanking region of the EAP1 gene was associated with increased incidence of amenorrhea (8,24). EAP1 has been shown to trans-activate the GnRH promoter but reduce Kiss1 transcription (25). This evidence together highlights the role of EAP1 in regulating GnRH expression and suggests that EAP1 may represent one of the main transcription factors contributing to the neuroendocrine control of female puberty.
Both of the mutations identified here in humans with DP led to an impaired ability of the mutant EAP1 protein to transactivate the human GnRH promoter, with a dose-dependent reduction in protein function seen in the in vitro heterozygous model (of 50% wild-type and 50% mutant protein). Interestingly, one of the mutations identified led to reduced protein function secondary to decreased expression of the EAP1 protein, while the second resulted in mislocalization of the mutant protein to the cytoplasm.
Disturbances in pubertal timing affect over 4% of adolescents and are associated with adverse health outcomes. Specifically, DP patients are at risk of decreased bone mineral density, osteoporosis (26,27) and psychological distress (28,29). In addition, recently associated increased risk has been shown for outcomes such as cervical cancer, myocardial infarction and poor overall health (12,30). Thus, the genetic control of human puberty is not only a fascinating scientific puzzle but also has potential for and transiently transfected with EAP1 plasmids (WT and mutated; 200 ng/well). Forty-eight hours post-transfection DLR assay was performed. Each transfection was normalized by co-transfecting with Renilla SV-40 vector and was performed in triplicate. EAP1 trans-activating GnRH promoter activity is significantly reduced by the in frame deletion (pAla221del) and missense (Asn770His) mutants compared to the WT. The mutants also cause a dose-dependent reduction of EAP1 WT transcriptional activity. The ablation of the RING finger domain (RINGdel) impairs the ability of trans-activating GnRH promoter and was used as a control, n = 3. (B) HEK293T cells were seeded at 0.3 × 10 6 cells/well onto a 6-well plate and transiently transfected with EAP1 plasmids (WT and mutated; 1 μg/well). Forty-eight hours post-transfection western blot analysis was performed. Eap1 protein expression is detected at the expected molecular weight (90 kDa) and RINGdel expression is detected at ∼80 kDa, as a result of the deletion of the RING finger domain (amino acids 715-762). GAPDH was used as loading control and detected at 37 kDa. HEK293T cells do not express EAP1, n = 3. (C) Quantification of western blot analysis, indicating that Ala221del protein levels are significantly reduced and Asn770His protein levels are significantly increased. major health impact for patients with abnormal pubertal timing, both precocious and delayed. Moreover, as the understanding of the genetic basis of both self-limited and other causes of DP improves, it is likely that genetic testing will be able to help to establish a definitive diagnosis in cases, particularly where there is diagnostic difficulty. Rapid and efficient diagnosis of patients in clinic would represent a significant advancement in patient care and a likely economic advantage.
In summary, we have identified two pathogenic mutations in the central transcriptional regulator EAP1 as the likely cause for self-limited DP in two families. Mutations in several genes known to influence GnRH secretion have been identified in patients with pubertal delay, including TAC3, TACR3, KISS1 and, now, EAP1 (31)(32)(33). Our results strengthen the concept that EAP1 is a bona fide regulator of pubertal onset and add to the understanding of the regulatory neural network that controls the onset of human puberty.

Patients
The cohort of individuals (n = 910) we studied here has been described in previous reports (14)(15)(16). Briefly, the cohort includes patients with self-limited DP (n = 492), defined as the onset of Tanner genital stage II (testicular volume, >3 ml) >13.5 years in boys or Tanner breast stage II >13.0 years in girls (i.e. 2 SD later than average pubertal development) (11) and their unaffected relatives. The patients were referred with a diagnosis of selflimited DP to a specialist pediatric care in Finland between 1982 and 2004. All patients met the diagnostic criteria for self-limited DP; medical history, clinical examination and routine labora-tory tests were reviewed to exclude those with chronic illness. Hypogonadotropic hypogonadism, if suspected, was excluded by spontaneous pubertal development at follow-up.
Family members of the DP patients participated via structured interviews and using archived height measurement records. The criteria for DP in probands' family members were (i) age at take-off or (ii) PHV occurring 1.5 SD beyond the mean (i.e. age at take-off exceeding 12.9 and 11.3 years or age at PHV exceeding 14.8 and 12.8 years in males and females, respectively) or (iii) age at attaining adult height more than 18 and 16 years, in males and females, respectively (34). Written informed consent was obtained from all participants. The study protocol was approved by the ethics committee for Pediatrics, Adolescent Medicine and Psychiatry, Hospital District of Helsinki and Uusimaa (570/E7/2003). UK ethical approval was granted by the London-Chelsea National Research Ethics Service (NRES) committee (13/LO/0257). The study was conducted in accordance with the guidelines of the Declaration of Helsinki.

Genetic analysis
Genetic analysis was performed in 67 probands with DP, from those 67 families with the greatest number of affected individuals in our cohort (male, n = 57; female, n = 10), 58 affected family members (male, n = 36; female, n = 22) and 35 of their unaffected family members (male, n = 13; female, n = 22). WES was performed on DNA extracted from peripheral blood leukocytes of these 160 individuals, using a Nimblegen V2 or Agilent V5 platform and Illumina HiSeq 2000 sequencing. The exome sequences were aligned to the UCSC hg19 reference genome using the Burrows-Wheeler Aligner software [BWA-MEM (bwa-0.7.12)]. The software Picard tools (picard-tools-1.119) was used to sort alignments and mark PCR duplicates. We used the genome analysis toolkit (GATK-3.4-46) to realign around indels and recalibrate quality scores using dbSNP, Mills and 1000 Genomes as reference resources. Variant calling and joint genotyping using pedigree information were performed using HaplotypeCaller in GVCF mode from the genome analysis toolkit. The resulting variants were filtered using the variant quality score recalibration from GATK.
Variants were analyzed and filtered for potential causal variants using filters for quality control, predicted function, MAF and biological relevance (Fig. 1A). Filtering by MAF included only variants with MAF <2.5% in the 1000 Genomes database, the NHLBI exome variant server and the Exome Aggregation Consortium (ExAC) and gnoMAD databases. Biological relevance filtering allowed prioritization of variants in genes with potential biological significance in the control of pubertal timing, using tools including Ingenuity Variant Analysis (QIAGEN Redwood City), Genego MetaCore (Thomson Reuters), OMIM, UniProt and Annovar (35). The multiple family filter retained only genes with variants present in more than one proband, and the case-control analysis excluded variants present in more than one unaffected control from our Finnish cohort, represented by family members with timing of puberty within the normal range. Targeted exome sequencing (Fluidigm) of the remaining candidate gene postfiltering was performed in a further 42 families from the same cohort [288 individuals, 178 with DP (male, 106; female, 72) and 110 controls (male, 55; female, 55); Fig. 1A], with filtering as in (15). Whole-gene RVBT was performed post-sequencing. Fisher's exact test was used to compare the prevalence of deleterious variants in our cohort with a set of controls from the Finnish population (ExAC European Finnish), taken from the ExAC Browser (ExAC, Cambridge, MA: accessed September 2015). For each gene, all variants from the ExAC database with a MAF of <2.5%, predicted to be deleterious by both Polyphen-2 (36) and SIFT (37), were included in the analysis, with each family in our cohort represented by the proband only. A multiple comparison adjustment was applied post-hoc using the Benjamini and Hochberg method (17), as detailed in (15).

In silico analysis
The mutant 3D structure of EAP1 harboring the Asn770His variant was obtained using Phyre2 program (38) and the mutant EAP1 FASTA sequence. Protein disorder was predicted using the DISOPRED3 server (39). The human experimental structure of EAP1 (PDB: 2cs3) was used as template.

Protein extraction and western blot analysis
Forty-eight hours post-transfection cells were harvested and lysed in radio-immunoprecipitation assay buffer (Sigma Aldrich) supplemented with protease inhibitor (Roche Diagnostics Ltd, West Sussex, UK). The concentration of the cell lysates was measured by the Pierce BCA protein assay kit (Thermo Fischer Scientific) according to the manufacturer's instructions. Equal amount of proteins was separated by standard deviation score (SDS)-PAGE (4-12% polyacrylamide NuPage Bis-Tris gels; Invitrogen) and transferred on nitrocellulose membranes (Promega, Southampton, UK). Non-specific binding was blocked with 5% non-fat milk in phosphate buffered saline (PBS) containing 0.1% Tween-20 (PBT). Membranes were incubated overnight at 4 • C with primary antibodies diluted in 5% non-fat milk in PBT: rabbit polyclonal anti-EAP1 diluted at 1:2000 (Sigma Genosys) and mouse monoclonal anti-GAPDH diluted at 1:5000 (Santa Cruz Biotechnology, Heidelberg, Germany). After washes in PBT, the membranes were probed with secondary antibodies (1:10 000; Licor, Cambridge, UK). After washes in PBT, the membranes were scanned and analyzed using the Odyssey Fc Imaging System (Licor). Experiment was repeated three independent times.

Promoter assay
HEK293T cells were plated on a 24-well plate (175 000 cells/well) and after 24 h were transiently co-transfected with the same amount of pGL4.10[luc]GnRH vector (7) and SV-40 Renilla (Promega, Southampton, UK) (150 ng/well each) in conjunction with either (i) EAP1 wild-type or mutant vectors (200 ng/well) or (ii) equal amount of EAP1 wild-type and p.Ala221del or wild-type and p.Asn770His mutants (100 ng/well each vector). The total amount of DNA transfected was kept constant to 500 ng/well by adding the appropriate amount of pBlueScript vector (i.e. in order to test the basal activation of the signal, co-transfection of pGL4.10[luc]GnRH and SV-40 Renilla only was also performed). Cells were harvested 48 h post-transfection and assayed for luciferase using Dual Luciferase Reporter System (Promega, Southampton, UK) following the manufacturer's instructions.
Each experiment was performed in triplicate and repeated three independent times. Samples were processed using POLARstar Omega microplate reader and data were analyzed using MARS Omega software.

In situ hybridization and immunohistochemistry
Peri-pubertal mouse brains were collected from timed crosses of C57BL/6 mice. Vaginal opening and balanopreputial separation were designated as indicative of the onset of puberty. Brains were fixed in 4% paraformaldehyde (PFA) in PBS, cryoprotected in 30% sucrose in PBS and frozen in OCT compound (VWR); 12 μm thick serial coronal sections were collected on Superfrost Plus slides (VWR). Mouse Eap1 was PCR-amplified from brain cDNAs using the following primers: Eap1 FOR: 5 -CAGTCTTGCTACCTGTGCGA-3 and Eap1 REV: 5 -AAGCGAGTGGTCCTTCTTGA-3 . Amplified cDNAs were cloned into the dual promoter vector pGEM-T easy (Promega, Southampton, UK) and linearized with the appropriate restriction enzymes. Probe preparation and in situ protocol were performed as previously in (40). When colabeling was desired, after in situ hybridization, the sections were incubated with primary antibodies (anti-GnRH; Immunostar, Hudson, Wisconsin, USA) diluted at 1:1000 in PBS-Triton 0.1%, as used in (41). After three washes with PBS-Triton 0.1%, slides were incubated with biotin-conjugated goat secondary antibodies (Vector Laboratories, Peterborough, UK), diluted at 1:300 in PBT and, after further washes, with the avidin-biotin complex (ABC staining kit; Vector Laboratories). The sections were reacted with 3,3 -diaminobenzidine (Vector Laboratories) and mounted in an aqueous compound formed by PBS and glycerol (3:1). Images were acquired using a Leica DM5500B microscope (Leica, Nussloch, Germany), equipped with a DCF295 camera (Leica) and DCViewer software (Leica), and then processed with Abode Photoshop CS6 and Adobe Illustrator CS6.

Non-human primates
The MBH of female rhesus monkeys (Macaca mulatta) was obtained through the Oregon National Primate Research Center Tissue Distribution Program. Animals were classified into different stages of pubertal development following the criteria proposed by Watanabe and Terasawa (42). EJ animals were aged from 9 months to 1.8 years and LJ animals were aged 2-2.9 years. The MBH was collected by making a cut along the posterior border of the optic chiasm, a cut in front of the mammillary bodies, and two lateral cuts half-way between the medial eminence and the hypothalamic sulci, as reported. Tissues were flash frozen in liquid nitrogen and stored at −80 • C for later processing.

ChIP assay
To study the recruitment of EAP1 to the monkey GnRH1 promoter, we performed ChIP assays using chromatin extracted from the MBH of two female rhesus monkeys. One of these animals (1 year and 65 days age) was in the early EJ phase, while the other (2 years and 99 days age) was in the LJ phase of postnatal development. The ChIP procedure was described previously (21,(43)(44)(45) and was carried out with slight modifications. Tissue cross-linking was performed by incubating the MBH fragments in 1% formaldehyde for 10 min at room temperature. After two additional washing steps in PBS, MBH fragments were lysed with 200 μl SDS buffer (0.5% SDS; 50 mm Tris-HCl; 10 mm EDTA) containing protease, phosphatase and HDAC inhibitors and sonicated for 45 s to yield chromatin fragments of ∼500 bp using a Fisher Scientific FB 705 sonicator. Size fragmentation was confirmed by agarose gel electrophoresis. The sonicated chromatin was clarified by centrifugation at 14 000 rpm for 10 min at 4 • C, brought up to 1 ml in ChIP dilution buffer (16.7 mm Tris-HCl, pH 8.1, 150 mm NaCl, 1.2 mm EDTA, 1.1% Triton X-100 and 0.01% SDS) containing protease, phosphatase and HDAC inhibitors and stored at −80 • C for subsequent immunoprecipitation. For this step, chromatin was pre-cleared with Protein A/G beads (Dynabeads; Invitrogen) for 1 h at 4 • C. One hundred microliter aliquots of chromatin were then incubated with 5 μg of the following antibodies: rabbit anti-EAP1 (custom made by Sigma) or anti-beta-galactosidase (Cortex, Biochem, Madison, WI). Antibody-chromatin complexes and 25 μl of protein A beads solution (Dynabeads; ThermoFisher, Waltham, MA) were incubated at 4 • C overnight with gentle agitation. Immunocomplexes were washed sequentially with 0.5 ml low-salt wash buffer (20 mm Tris-HCl, pH 8.1, 150 mm NaCl, 2 mm EDTA, 1% Triton X-100 and 0.1% SDS), high-salt wash buffer (20 mm Tris-HCl, pH 8.1, 500 mm NaCl, 2 mm EDTA, 1% Triton X-100 and 0.1% SDS), LiCl buffer (10 mm Tris-HCl, pH 8.1, 250 M LiCl, 1% Nonidet P-40, 1% sodium deoxycholate and 1 mm EDTA) and TE buffer (10 mm Tris-HCl, pH 8.0 and 1 mm EDTA). The immunocomplexes were eluted with 100 μl of 0.1 M NaHCO3 and 1% SDS at 65 • C for 45 min. Cross-linking was reversed by adding 4 μl of 5 M NaCl and incubating at 95 • C for 30 min. DNA was recovered by using the ChIP DNA clean and concentrator columns (Zymo Research, Irvine, CA) and stored at −80 • C until subsequent PCR analysis. All chemicals were purchased from Sigma-Aldrich. Due to the limitation of availability of primate tissue, the experiment was performed only once.

PCR detection of chromatin-immunoprecipitated DNA
The promoter region of the rhesus monkey GnRH gene (JU472587) was amplified using the forward primer GCCA-GAAGCTTCCAGACATCC and the reverse primer AAGTGCAGC-CATTAAAACCTCAG. As a negative control for EAP1 binding we used an intergenic region located in intron 2 of the GnRH1 gene (forward primer: ACCACGCCCGGACTGTTTC and reverse primer: TGATCCACTTACCTCGGCTTCC; Eurofins MWG Operon, Huntsville). PCR reactions were performed using 1 μl of each IP and input samples and HotStart Taq Polymerase (Qiagen, Germantown, MD) in a volume of 25 μl. The thermocycling conditions used were 95 • C for 5 min, followed by 33 cycles of 15 s at 95 • C followed by 30 s at 55 • C and 30 s at 72 • C. The PCR products were run in a 1.2% agarose gel prepared in Tris/Borate/EDTA buffer.

Immunofluorescence
Cells were seeded in a 4-well cell culture slide at a density of 1 × 10 5 cells/well (Millipore, Fisher Scientific, Watford, UK) and transiently transfected with 500 ng/well of WT or mutant Ala221del and Asn770His expression plasmids using PEI. Fortyeight hours after transfection, cells were fixed in 4% PFA in PBS for 15 min and washed with PBS. Samples were permeabilized with 0.1% Triton X-100 in PBS for 30 min and blocked with blocking buffer (10% normal goat serum in PBS) for 30 min. The staining was performed by incubating the samples with rabbit anti-EAP1 diluted at 1:1000 (Sigma Genosys) in blocking buffer for 1 h, followed by a 30 min incubation with goat anti-rabbit Alexa fluor 594 (ThermoFisher Scientific; 1:250) antibody. The cell nuclei were stained with DAPI (Sigma). Images were acquired using a fluorescence microscope (Leica microsystem, Germany) and processed using Adobe Photoshop CS6.

Statistical analysis
For all experiments, data are expressed as the mean ± SEM. One-way repeated measures analysis of variance was used to determine statistical significance for multiple comparisons; Pvalues ( * ) of <0.05 and ( * * ) of <0.01 were considered statistically significant. A P-value ( * * * ) of <0.001 was considered highly significant. The statistical analysis was performed using GraphPad Prism7 (GraphPad software).