An Association between EMX2 Variations and Mayer-Rokitansky-Küster-Hauser Syndrome: A Case-Control Study of Chinese Women.

Mayer-Rokitansky-Küster-Hauser (MRKH) syndrome is characterized by congenital malformations of Müllerian structures, including the uterus and upper two-thirds of the vagina in women. Until now, the etiology of this disease has remained unknown. We hypothesized that EMX2 (the human homologue of Drosophila empty spiracles gene (2) might be a candidate gene for MRKH syndrome because it plays an important role in the development of the urogenital system. Through sequence analysis of EMX2 in forty patients with MRKH syndrome and one hundred and forty healthy women controls, we identified eleven variations in total. Four novel variations were only found in MRKH patients, and seven single nucleotide polymorphisms were identified in both patients and controls. In silico analyses suggested that the novel variations in the 5'UTR (untranslated region) and 3'UTR might affect transcriptional activity of the EMX2 promoter or posttranscriptional processing. In conclusion, our study suggests an association between noncoding variations in the EMX2 gene and MRKH syndrome in a Chinese Han population.

1. Introduction

Mayer-Rokitansky-Küster-Hauser (MRKH) syndrome is characterized by congenital malformations of Müllerian structures, including the uterus and upper two-thirds of the vagina. It is also referred to as CAUV (congenital absence of the uterus and vagina) or MA (Müllerian aplasia). This syndrome is a rare disease, affecting 1 in every 4500 female live births [1]. MRKH patients usually have a 46, XX karyotype, normal secondary sexual characteristics and normal ovaries with no sign of androgen excess. Primary amenorrhea is the most prominent feature in these patients [2]. Other malformations that are often associated with this disorder include renal malformations, skeletal abnormalities, hearing defects, and heart malformations [3].

Although the majority of MRKH cases are sporadic, familial occurrence suggests a genetic cause. Thus, a candidate gene approach has been adopted based on genes involved in Müllerian duct development during embryogenesis or other genetic diseases that share characteristics similar to MRKH syndrome. These candidate genes, such as WNT family members [4-6], HOXA family members [7, 8], TCF2 [9], PAX2 [10], and LHX1 [6], have been studied in patients with MRKH syndrome. Unfortunately, these results are unproductive, and the molecular basis of MRKH syndrome remains to be elucidated. Although WNT4 mutations were identified in certain patients, all of these patients presented with uterovaginal aplasia and hyperandrogenism, which is considered to be distinct from MRKH syndrome [6, 11, 12]. Moreover, heterozygous mutations of LHX1 were detected in two patients, and it was suggested that LHX1 mutations might be the cause of MRKH syndrome in a subgroup of patients [13]. Additionally, it has been widely accepted to date that MRKH syndrome might occur due to polygenic/multifactorial inheritance [3, 14].

Molecular expression analyses of targeted mutagenesis in mouse models have helped to identify several genes that are involved in the development of the urogenital system. Emx2 (empty spiracles homeobox (2) is a divergent homeobox-containing gene orthologous to the Drosophila empty spiracles gene (ems) and is involved in the development of the mammalian brain and urogenital system. Mammalian embryos have both Wolffian ducts and Müllerian ducts, and Müllerian ducts normally develop in parallel with Wolffian ducts approximately 13.0 dpc (day after coitum) in both male and female wild-type mice [15]. Müllerian ducts then differentiate into the oviducts, uterus, cervix, and upper portion of the vagina in females, while they degenerate in males. Emx2 is expressed in the epithelial cells of Wolffian and Müllerian ducts [16]. In homozygous Emx2 mutant mice, the kidneys, ureters, gonads, and genital tracts are completely absent, and Müllerian ducts never form. The phenotype of Emx2 mutant mice is similar to that of MRKH patients. We therefore postulated a connection between the EMX2 gene and MRKH syndrome. We performed a sequence analysis of EMX2 variants in a case-control study. To the best of our knowledge, our study was the first to examine whether MRKH syndrome occurs due to variations in the EMX2 gene in hope of elucidating the pathogenesis of MRKH syndrome.

2. Materials and Methods

2.1. Subjects

Forty Han Chinese patients with sporadic MRKH syndrome and one hundred and forty randomized matched female controls with a normal reproductive history, i.e., at least one normal pregnancy without history of genital abnormalities, were analyzed. All patients were included in this study according to the following criteria: 46, XX karyotype, normal secondary sexual characteristics, primary amenorrhea, and absence of uterus, cervix, and proximal vagina documented by ultrasonography and laparoscopy. Six patients presented with renal abnormalities, including unilateral renal aplasia (five patients) and renal ectopia (one patient). Of the patients with renal abnormalities, one had inguinal hernia and another had thoracic vertebral malformations. This study was approved by the local ethics committee, and informed consent was obtained from each participant before entry into the study.

Total genomic DNA from all participants was isolated from peripheral blood using a QuickGene DNA whole blood kit S (Fujifilm, Japan).

2.2. Polymerase Chain Reaction (PCR) and Sequencing

All three exons and exon-intron boundaries of the human EMX2 gene [reference sequence NM_ 004098.3, hg 19/GRCh 37] were amplified by PCR. The primers were designed using Primer Premier 5 software and are presented in Table 1. Each PCR reaction was performed in a total volume of 25 μl containing 20–50 ng of genomic DNA, 0.4 μM each primer, 10 × LA PCR buffer II (Mg2+ plus), 0.4 mM dNTP mixture, and 1.25U of LA Taq DNA polymerase (TaKaRa, Japan). A touchdown PCR program was used to amplify exon 1, whereas conventional PCR reactions were used for exon 2 and exon 3. The PCR conditions are shown in Table 1. PCR products were subsequently sequenced and analyzed using an ABI 3730XL Genetic Analyzer (Applied Biosystems, USA). To verify whether the detected variations were located in conserved regions, Genome Browser of UCSC (http://genome.ucsc.edu/) was used to render sequence alignments in different species, and constraint scores were calculated by genomic evolutionary rate profiling (GERP) [17, 18]. The DataBase of Transcriptional Start Sites (DBTSS, http://dbtss.hgc.jp) and the MatInspector program were used for transcriptional start site analysis and the identification of binding sites for transcription factor prediction, respectively [19]. Moreover, we determined the effect of genetic variations on RNA folding using GeneQuest and analyzed possible microRNA binding sites using TargetScanHuman 6.2 (http://www.targetscan.org).

Table 1

Primers and PCR conditions for EMX2 gene amplification.

Primer	Primer sequence (5′-3′)	PCR conditions
Exon1F	CGCTAGGCTAGAGGAATCTGTCTGT	3 min at 97°C; 30 s at 97°C, 30s at 68°C (−0.5 °C per cycle), 90 s at 72°C (30 cycles), and 3 min at 97°C; ‘30 s at 97°C, 30 s at 61°C, 90 s at 72°C (5 cycles), and 10 min at 72°C
Exon1R	CAACTTCTCCGTTCGCACCC
Exon2F	GAGCAGGCGTTCCCTTCGT	3 min at 94°C; 1 min at 94°C, 1 min at 63.4°C, 1 min at 72°C (35 cycles), and 10 min at 72°C
Exon2R	GGGCAGATTGGCACTTACAGC
Exon3F	GCTGGGTCTTTGCTGAGTC	3 min at 94°C; 1 min at 94°C, 45s at 67.3°C, 2 min at 72°C (35 cycles), and 10 min at 72°C
Exon3R	TGAGGAGCCTGGGTTTCTT

2.3. Statistical Analysis

All statistical analyses were conducted using the SPSS statistical package (version 17.0). The Hardy-Weinberg equilibrium for each of the variations in controls was assessed based on a goodness-of-fit χ2 test. The differences of allelic and genotypic distributions between patients and controls were measured using a χ2 test. Logistic regression analysis was used to calculate the odds ratio and 95% confidence interval (95% CI) values. The level of significance was taken as P < 0.05. Primers and PCR conditions for EMX2 gene amplification ARE shown in Table 1.

3. Results

By sequencing the entire coding region, exon-intron boundaries, 5′UTR and 3′UTR of EMX2, we identified a total of eleven variations; the results are shown in Tables 2 and 3.

Table 2

EMX2 variations detected only in MRKH patients.

Exon/intron	Position	Genotype	No. (%) in patients
Exon 1	c.252 A > G	AG	1 (2.5)
		AA	39 (97.5)
Exon 1	c.-433_-432insC	insC/wt	1 (2.5)
		Wild/wt	39 (97.5)
Exon 1	c.-621G > C	GC	1 (2.5)
		GG	39 (97.5)
Exon 3	c.∗950C > T	CT	1 (2.5)
		CC	39 (97.5)

Table 3

EMX2 variations detected both in MRKH patients and controls.

Exon/intron	dbSNP	Position	GERP score	Genotype	No. (%) in patients	No. (%) in controls
Exon 1	rs12777466	c.-674G > A†	5.07	AA	2 (5)	14 (1)
				AG	16 (40)	37 (26.4)
				GG	22 (55)	89 (72.6)
Intron 1	rs8192644	c.406 + 113G > A	2.04	AA	0 (0)	1 (0.7)
				GA	4 (10)	9 (7.5)
				GG	36 (90)	130 (91.8)
Intron 1	rs142080828	c.407–104G > A	1.23	GA	1 (2.5)	3 (2.1)
				GG	39 (97.5)	137 (97.9)
Intron 2	rs202171958	c.591+6T > C	−6.74	TC	1 (2.5)	1 (0.7)
				TT	39 (97.5)	139 (99.3)
Exon 3z	rs66710107	c.∗317_∗320del	—	del/del	40 (100)	40 (100)
		AGAG		wt/wt	0 (0)
Exon 3	rs187010704	c.∗925C > T	4.88	CT	1 (2.5)	1 (0.7)
				CC	39 (97.5)	139 (99.3)
Exon 3	rs41284394	c.∗1201T > C	5.79	CC	2 (5)	3 (2.1)
				CT	0 (0)	5 (3.6)
				TT	38 (95)	132 (94.3)

Four of these variations were only detected in patients, and the others were present both in patients and controls. The allele and genotype frequencies for all the variants were determined based on the Hardy-Weinberg equilibrium.

We identified four heterozygous variations (c.-621G > C, c.-433_-432insC, c.252 A > G, and c.950C > T) that were only present in MRKH patients, each with an allelic frequency of 1.25%; the results are shown in Table 2. All of these variations were novel and were not annotated in dbSNP141. The first two variations (c.-621G > C and c.-433_-432insC) were located in the 5′UTR of exon 1. The sequence change c.-621G > C (GERP score: 5.04) is located in a highly conserved region, as verified using the UCSC Genome Browser tool. According to DBTSS, the position of this variation was a possible transcriptional start site of the EMX2 gene. Analysis using the MatInspector program revealed that the region was possibly associated with several transcriptional factors, such as homeobox transcription factors. The second variation, c.-433_-432insC, was not located in a highly conserved region, although the region was predicted to interact with at least one of the three specific transcription factors: HMX2, MSX1, and MSX2. The synonymous nucleotide substitution, c.252 A > G, was also not located in a conserved region. Another change, c.∗950C > T, is located in the 3′UTR of exon 3, and analyses using the UCSC Genome Browser tool and GERP test (GERP score: 3.94) suggested that the locus of this variation was conserved among mammals. Analysis using GeneQuest showed that this variation would result in distinctive RNA folding compared with the reference sequence.

Moreover, seven previously reported polymorphisms were found both in patients and controls: one deletion and six single nucleotide polymorphisms, none of which were in coding regions; the results are shown in Table 3. However, no statistically significant associations were found for theses polymorphisms at the allele or genotype levels. We identified six single nucleotide polymorphisms: rs12777466 in exon 1, rs8192644 and rs142080828 in intron 1, rs202171958 in intron 2, and rs187010704 and rs41284394 in exon 3. The first variation, rs12777466, is located in a highly conserved region, and its GERP score of 5.07 was an indication of evolutionary conservation. However, the intronic variations rs8192644, rs142080828, and rs202171958 were not within highly conserved sequences and had GERP scores of 2.04, 1.23 and, −6.74, respectively.

The UCSC Genome Browser tool showed that the position of rs187010704 (c. ∗ 925C > T) was conserved in other mammals, except for mouse and rat and that of rs41284394 (c.∗1201T > C) was highly conserved in various species, such as zebrafish and Xenopus tropicalis. Both of these variations had high GERP scores of 4.88 and 5.79, respectively. Additionally, TargetScanHuman 6.2 analysis predicted that rs41284394 was within the possible binding site for miR-181abcd/4262. Additionally, the variation might change the RNA-folding process, according to the GeneQuest results. However, rs187010704, the other variation in the 3′UTR, was neither within any possible binding sites for microRNAs nor had an effect on RNA folding. In addition, we detected one deletion, rs66710107 in exon 3, that was not located in a conserved region based on the analysis of multiple sequence alignments.

4. Discussion

Patients with MRKH syndrome suffer from infertility and psychological distress. The disease results in high costs not only to the patient herself but also the whole family. Unfortunately, the etiology, prenatal diagnosis, and effective treatment of MRKH syndrome are currently unavailable. To identify genetic risk factors in MRKH patients, potential pathogenic mechanisms and susceptibility to the disease must be investigated.

There have been many unsuccessful attempts to identify genetic risk factors, and to the best of our knowledge, the implication of EMX2 variations in MRKH syndrome has not been studied. Emx2 plays an important role in urogenital development [15, 20]. In addition, Emx2 mutant mice do not have Müllerian ducts and present a phenotype similar to that observed in MRKH patients. Moreover, Emx2 is also crucial in the morphogenesis of the central nervous system and inner ear development [21, 22]. There are reports of MRKH patients with learning disabilities, mental impairment, hearing loss, and endometriosis, reinforcing the link between the EMX2 gene and MRKH syndrome [3, 23, 24]. To investigate the possibility of a link between EMX2 variations and MRKH syndrome, we screened forty patients with MRKH syndrome and one hundred and forty healthy females for variations in the EMX2 gene. A total of eleven variations were identified.

We detected four heterozygous variations (c.-621G > C, c.-433_-432insC, c.252 A > G, and c.∗950C > T) in MRKH patients. The allele frequency of these variations was 1.25% in MRKH patients, and none of these variations were found in our matching control cohort or dbSNP141. It should be noted that all of these variations are reported in our study for the first time. The 5′UTR variations c.-621G > C and c.-433_-432insC might affect the transcriptional regulation of EMX2. According to in silico analysis using MatInspector, c.-621G > C is located within the potential binding sites for transcriptional factors (including HOXA1, HOXA5, HOXA7, HOXB1, HOXB4, HOXB5, HOXC6, HOXD3, HOXD4, EN2, GSX2, and VAX1), and c.-433_-432insC might be associated with transcription factors of HMX2, HOX7, or HOX8. Interestingly, EMX2 is negatively regulated by HOXA10 binding to a 150 bp EMX2 regulatory element [25, 26]. Furthermore, HOX genes are known to encode transcription factors that play crucial roles in the development of the female reproductive tract. These two 5′UTR variations are not located within the known binding site for HOXA10, although their locations were predicted to interact with other transcriptional modulators, including HOX genes. The HMX2 gene plays an important role in organ development during embryogenesis, especially in inner ear formation [27, 28]. Moreover, bioinformatic analysis showed that c.-621G > C is located in a highly conserved region.

Based on the conservation analyses, this region might be of great importance to the function of the HMX2 gene. It is also possible that the c.-621G > C variation might change the transcriptional start site of the EMX2 gene based on the DBTSS analysis. These variations may change the promoter function and therefore affect the expression of the EMX2 protein. Another synonymous variation, c.252 A > G, is not located in a conserved region and is therefore unlikely to be pathologic. Currently, the important role of the 3′UTR in posttranscriptional regulation, such as mRNA stability and degradation, is well understood [29]. The 3′UTR variation in the EMX2 might also have an effect on gene function. The variation c.∗950C > T is located in a highly conserved region according to the bioinformatic analysis. When compared with the reference sequence, c.∗950C > T resulted in a distinct RNA-folding process. Overall, these variations in regulatory domains might alter the transcriptional activity of the EMX2 promoter or have posttranscriptional effects. Although the low frequency of these variations may be due to small sample size, it should be validated in larger MRKH patient groups.

We also detected seven other previously reported polymorphisms (rs12777466, rs8192644, rs142080828, rs202171958, rs66710107, rs187010704, and rs41284394). The variation rs12777466 is located in the 5′UTR, 674 bp before the start codon in the EMX2 gene. UCSC Genome Browser and GERP score analyses showed that this variation occurs in a conserved region. Its location might be associated with another transcription factor, but not the known binding site for HOXA10. We could not exclude the possibility that the associated transcription factor plays a crucial role during reproductive tract development. Variations in the 5′UTR might affect gene expression at the translational level, and a high level of conversation indicates that this nucleotide residue might be quite important. Three variations (rs8192644, rs142080828, and rs202171958) in the introns and one deletion rs66710107 (c.∗317_320delAGAG) in the 3′UTR were not found to be located in conserved regions. However, the bioinformatic analysis showed that two variations in the 3′UTR of exon 3, rs187010704, and rs41284394, are located in highly conserved regions. It was predicted that rs41284394 lies within the possible binding site for miR-181abcd/4262. Interestingly, previous studies suggested important roles of miR-181 in osteoblastic differentiation and the immune system [30, 31]. Additionally, skeletal abnormalities, such as scoliosis and vertebral anomalies, are associated with MRKH syndrome [1]. Thus, miR-181 might be a link between EMX2 and MRKH syndrome. The GeneQuest analysis also suggested that rs41284394 might alter RNA folding, which might be involved in posttranscriptional regulation. However, the allele and genotype frequencies of all of these variations were not significantly different between MRKH patients and controls. This might be due to our small sample size, and further experimental investigations are needed to determine whether these variations are linked to MRKH syndrome.

5. Conclusions

In conclusion, four variations were reported for the first time in our study, all of which were absent in our control group and dbSNP141. According to in silico analyses, two of these variations (c.-621G > C and c.∗950C > T) in the EMX2 gene might be associated with an increased susceptibility to MRKH syndrome in a Chinese Han population. This study provides the first insight into the involvement of the EMX2 gene in genetic predisposition to MRKH syndrome. Our results do not suggest associations between coding variations in the EMX2 gene and MRKH syndrome. However, we cannot exclude the possibility that regulatory variations in the 5′UTR and 3′UTR might contribute to the polygenic/multifactorial pathogenesis of MRKH syndrome. These findings should be confirmed in large-scale studies, and in vitro functional studies are needed for further evaluation of the association of EMX2 variations with MRKH syndrome.