NRG1 variant effects in patients with Hirschsprung disease

Background Hirschsprung disease (HSCR) is a heterogeneous genetic disorder characterized by absence of ganglion cells along the intestines resulting in functional bowel obstruction. Mutations in neuregulin 1 (NRG1) gene have been implicated in some cases of intestinal aganglionosis. This study aims to investigate the contribution of the NRG1 gene to HSCR development in an Indonesian population. Methods We analyzed the entire coding region of the NRG1 gene in 54 histopathologically diagnosed HSCR patients. Results All patients were sporadic non-syndromic HSCR with 53/54 (98%) short-segment and 1/54 (2%) long-segment patients. NRG1 gene analysis identified one rare variant, c.397G > C (p.V133 L), and three common variants, rs7834206, rs3735774, and rs75155858. The p.V133 L variant was predicted to reside within a region of high mammalian conservation, overlapping with the promoter and enhancer histone marks of relevant tissues such as digestive and smooth muscle tissues and potentially altering the AP-4_2, BDP1_disc3, Egr-1_known1, Egr-1_known4, HEN1_2 transcription factor binding motifs. This p.V133 L variant was absent in 92 non-HSCR controls. Furthermore, the rs7834206 polymorphism was associated with HSCR by case–control analysis (p = 0.037). Conclusions This study is the first report of a NRG1 rare variant associated with HSCR patients of South-East Asian ancestry and provides further insights into the contribution of NRG1 in the molecular genetic pathogenesis of HSCR. Electronic supplementary material The online version of this article (10.1186/s12887-018-1265-x) contains supplementary material, which is available to authorized users.


Background
Hirschsprung disease (HSCR), a heterogeneous genetic disorder, is characterized by the lack of ganglion cells along varying lengths of the intestines resulting in functional obstruction among children [1]. According to the length of aganglionosis, HSCR is categorized into three major types: short-segment (aganglionosis up to the upper sigmoid colon), long-segment (aganglionosis up to the splenic flexure and beyond) and total colonic aganglionosis (TCA) [1,2].
The incidence of HSCR varies among populations with 15, 21 and 28 cases per 100,000 live births in Caucasians, Africans and Asians, respectively [1,2]. These differences might be influenced by susceptibility factors such as the REarranged during Transfection (RET) rs2435357 risk allele frequency [3]. Our recent studies clearly demonstrated the presence of a higher frequency of RET rs2435357 susceptibility allele in Indonesian ancestry (50%) [4,5]. Nevertheless, we observed that semaphorin 3 rs11766001 common variant had different effects on HSCR depending on the ethnic background [6].
HSCR is a complex genetic disorder and at least 15 genes have been implicated in the development of HSCR, with the gene encoding the receptor tyrosine kinase RET accounting for up to 21% of sporadic cases [7]. These genes encode proteins that are important for the enteric ganglia development and are classified into three groups: 1) those that are associated with RET pathways (RET, GFRα1, GDNF, NTN, PSPN); 2) those implicated in endothelin type B receptor/EDNRB (EDNRB, EDN3, ECE-1) pathways; and 3) transcription factors that influence both RET and/or EDNRB pathways (SOX10, ZFXH1B, PHOX2B) [7,8].
Neuregulin 1 (NRG1) has been implicated as a disease susceptibility gene in Chinese [9,10]. Both common (single nucleotide polymorphisms, SNPs) and rare variants of this gene are reported to confer disease risk with over-representation in patient cohorts. However, in the first study on Caucasian, Luzón-Toro et al. [11] did not find any significant association of NRG1 SNPs in Spanish HSCR suggesting that the association of such common polymorphisms to the disease may be restricted to specific populations. Furthermore, they found three novel NRG1 rare variants that were hypothesized to have functional consequences during embryonic development of HSCR which could lead to HSCR in their patients. It is believed that NRG1 plays an essential role in the signaling pathway of the enteric nervous system (ENS) [9][10][11][12]. In addition, the development of ENS requires a balance of neurogenesis and gliogenesis by RET/GDNF and ERBB2/NRG1 pathways, respectively. Nrg1 suppresses Gdnf-induced neuronal differentiation and Gdnf negatively controls Nrg1-signaling by reducing the expression of its receptor, ErbB2 [13]. Recently, we studied two previously reported NRG1 SNPs and demonstrated that that only one common variant rs7835688 as a genetic risk factor for Indonesian HSCR [4]. Therefore, we aimed to perform a comprehensive study of the entire NRG1 gene to investigate fully if there is any link between NRG1 variants with Indonesian HSCR patients.
The study was reviewed and approved by the Institutional Review Board (IRB) of the Faculty of Medicine, Universitas Gadjah Mada/Dr. Sardjito Hospital, Indonesia (KE/FK/787/ EC/2015). Written informed consent was obtained from all parents of the subjects for this study.
Genomic DNA extraction and polymerase chain reaction (PCR) Genomic DNA was extracted from whole blood from each individual by using QIAamp DNA Extraction Kit (QIAGEN, Hilden, Germany), according to the manufacturer's instructions. The extracted DNA samples were stored at − 20°C until analysis. PCR was conducted by using a Swift Maxi thermal cycler (Esco Micro Pte. Ltd., Singapore). The primer sequences for NRG1 gene analysis were designed according to a previous study [10].

DNA sequencing
To identify mutations, all 16 exons of the NRG1 that had been PCR amplified in 20 fragments were screened using Sanger sequencing analysis with a BigDye Terminator V3.1 Cycle Sequencing Kit (Applied Biosystems, Foster City, CA). The products were separated and analyzed on a 3730xl Genetic Analyzer (Applied Biosystems, Foster City, CA) using DNA Sequencing Analysis Software (Applied Biosystems, Foster City, CA) [17,18].
SIFT predicts whether a change of amino acid change in a protein will have an impact on phenotype. This calculation is based on the premise that protein evolution is correlated with protein function, the more important positions of amino acids will correspondingly be more conserved in an alignment of the protein family. A SIFT score ranges from 0 to 1 and is classified as two groups: 1) predicted damaging if score ≤ 0.05; and 2) tolerated if score > 0.05 (http://sift.jcvi.org/). PolyPhen-2 score calculates the possible effect of a change of amino acid on the human protein structure and function. Its score ranges from 0 to 1.0 and consists of three interpretation: 1) benign (≤ 0.452 in PolyPhen2 HDiv and ≤ 0.446 in PolyPhen2 HVar); 2) possibly damaging (0.453 to 0.956 in PolyPhen2 HDiv and 0.447 to 0.909 in PolyPhen2 HVar); and 3) probably damaging (≥ 0.957 in PolyPhen2 HDiv and ≥ 0.909 in PolyPhen2 HVar) (http://genetics. bwh.harvard.edu/pph2/). Mutation Taster determines the disease potential of an alteration of amino acid using the Bayes classifier, with four possible analysis: 1) disease causing (i.e. probably deleterious); 2) disease causing automatic (i.e. known to be deleterious); 3) polymorphism (i.e. probably harmless); and 4) polymorphism automatic (i.e. known to be harmless) (http://www.mutationtaster.org). DANN uses deep neural network (DNN) algorithm to predict pathogenicity of variants and scores from 0.98 to 1 are considered to be protein disrupting (Additional file 1). All these tools predict whether a mutation leads to damaging function that may be potentially disease causing but it is well acknowledged that such in silico tools are only used as predictive guides and may not necessarily agree with each other as they can yield conflicting results [19][20][21].

Statistical genetic analysis
Chi-square test was used to establish p-values for the case-control association analysis for NRG1 variants (rs7834206, rs3735774, and rs75155858). A p-value less than 0.05 was considered significant.
We compared the observed NRG1 rare variant allele frequency in Indonesian HSCR patients with those reported for the 1000 Genomes Project and ExAC ancestry controls. The frequency of c.397G > C (p.V133 L) rare variant in our HSCR patients (0.9%) was higher than those reported for the 1000 Genomes and ExAC    ancestry controls (0 and 0.07%, respectively, p < 0.0001) ( Table 3). Next, we analyzed the potential damaging effect of the rare variant (p.V133 L) using the following function prediction algorithms: SIFT, PolyPhen-2 (HDiv and HVar), LRT, Mutation Taster, Mutation Assessor, FATHMM, CADD and DANN (Table 4). Although SIFT and PolyPhen-2 analysis of p.V133 L showed the variant as being tolerated (0.22) and benign (0.029 in PolyPhen2 HDiv and 0.02 in PolyPhen2 HVar), respectively, p.V133 L was predicted to be disease causing by Mutation Taster and DANN. According to the conservation scores prediction using GERP and PhyloP vertebrate, the p.V133 L variant reached a deleterious threshold with a score of 4.26 and 1.799, respectively (Table 5) [22,23]. However, ClinVar reported the clinical significance of the p.V133 L variant as likely benign (Table 6).
We then used the HaploReg database to assess the regulatory potential of this variant. This variant was predicted to reside within a region of high mammalian conservation, overlapping with the promoter and enhancer histone marks of relevant tissues such as digestive and smooth muscle tissues and altering AP-4_2, BDP1_disc3, Egr-1_known1, Egr-1_known4, HEN1_2 transcription factor binding motifs (Table 7). When 92 Indonesian control subjects were screened, none carried this variant, suggesting that it is not likely to be a common variant in our population.

Discussion
In this study, we have conducted a mutational screening of the NRG1 gene in an Indonesian cohort of HSCR patients. We found three non-synonymous variants and one nucleotide substitution variant in the 5′ untranslated region of the NRG1 gene. Among them, p.V133 L located in the EGF-like domain within the NRG1 protein, was considered a rare variant in our population since it was absent in 92 non-HSCR controls. There was conflicting results from computational prediction programs with regards to the pathogenicity of p.V133 L and it has been classified as benign in ClinVar. However, the location of this variant in a well conserved and important histone and motif binding region suggests that it may have a regulatory function on gene expression.
The EGF domain where p.V133 L resides has been shown to be necessary for the activation of ErbB receptors [24]. The variant might also affect the SOX10-mediated maintenance of ENS progenitors since the receptors for Nrg1 and ErbB3 are regulated by Sox10 [25,26]. Moreover, the ENS development involves a balance of neurogenesis and gliogenesis by RET/GDNF and ERBB2/NRG1 pathways, respectively [13]. Therefore, we might suggest that p.V133 L will have an effect in intestinal aganglionosis in our HSCR patient by two mechanisms: 1) affecting the SOX10-mediated maintenance of ENS progenitors, and 2) altering the balance of neurogenesis and gliogenesis during ENS development. The p.V133 L variant was absent in 92 non-HSCR population. Thus, we considered the variant as rare in our population. The presence of this variant was also higher in our patients than those reported in the 1000 Genomes and ExAC ancestry controls (0.9% vs. 0 vs. 0.07%). Previous studies reported NRG1 rare variants to be implicated in the development of ENS and HSCR [10,11]. Interestingly, the frequency of this NRG1 rare variant found in Indonesian [1/54 (1.8%) patients] was similar to those of Caucasian and Chinese probands [3/207 (1.4%) and 13/358 (3.6%) patients, respectively; with p = 0.45] [10,11]. Furthermore, the NRG1 rs7835688 common variant was originally associated with HSCR in Chinese patients [9] and has been similarly observed in other Asian ancestry cases [4,27], but these associations were not replicated in any Caucasian population [11,28].
We were able to find evidence of the genetic effect of NRG1 rs7834206 in Indonesian HSCR cases. Our study shows that NRG1 rs7834206 is a genetic risk factor for HSCR with a background allele frequency of~19% in Indonesian populations ( Table 3). The risk allele at rs7834206 (0.19 vs. 0.18) in Indonesian controls has a similar frequency to subjects of Asian ancestry in the 1000 Genomes dataset (p = 0.80). Furthermore, a recent meta-analysis study showed that NRG1 polymorphisms are risk factors for HSCR in Asians but not in Caucasians [29]. Therefore, our findings strengthen the notion of the potential damaging role of NRG1 common variants in HSCR patients of Asian ancestries. It should be noted that the small sample size in this study poses a potential limitation and a significantly larger number of patients are needed to validate these observations. Future investigations on larger cohorts will clarify the exact role of the two variants, rs7834206 and p.V133 L in HSCR in Indonesian population. Moreover, moving beyond the prediction effect of mutations, further in vitro or in vivo functional studies are required to shed light on the actual effect of any variant.
Although there have been previous studies correlating NRG1 SNPs with disease association in HSCR, these have only been reported in three populations, namely Chinese, Caucasian, and Thais [10,11,27,28]. Indeed, contrary to the findings in Chinese, Luzón-Toro et al. [11] failed to find association of NRG1 variants with HSCR phenotype in Spanish patients and suggested that this could be due to population differences. In the Thai patients, only four NRG1 SNPs were selected for association studies and disease correlation was found in Thai-Chinese but not Thai-Muslim patients suggesting ethnicity differences [27]. Our present study is unique in this being a comprehensive screen of NRG1 gene in Indonesian HSCR patients. As indicated in all the above studies, it is important to verify the role of NRG1 variants in different population groups due to ethnicity differences. Furthermore, the SNPs in this study, namely the three common variants (rs7834206, rs3735774, and rs75155858) have not been previously investigated in all the other reported studies [10][11][12]. Since the allele frequencies of some common variants have been known to vary among different ethnic groups within Asian population, it is thus important to clarify the relationship of NGR1 SNPs in Indonesian patients [30].
In view of recent interest in the role of NRG1 in HSCR, it is not surprising that there is increasing evidence showing the contributions of both common and rare variants of this gene to disease risk. Our study reinforces this growing body of information showing the presence of NRG1 risk variants that may interact with other alleles such as RET, GDNF, ErbB or SOX10 as well as other susceptibility factors in leading to HSCR risk, and demonstrating that there are population differences in the contributions of the different risk alleles.

Conclusion
This study is the first report of a NRG1 rare variant associated with HSCR patients in South-East Asian ancestry and adds insights into the role of NRG1 in the molecular pathogenesis of HSCR.