Association between IRF6 and 8q24 polymorphisms and nonsyndromic cleft lip with or without cleft palate: Systematic review and meta‐analysis

Background We conducted a systematic review and meta‐analysis of interferon regulatory factor 6 and 8q24 polymorphisms with nonsyndromic cleft lip with/without cleft palate (NSCL/P). Methods Data extraction was independently performed by two reviewers. Genotypic effects of four polymorphisms from 31 studies were pooled separately by ethnicity using a mixed‐effect logit model with accounting for heterogeneity. Results For rs2235371, AA and GA carried, respectively, 51% (95% confidence interval [CI], 37%–61%) and 42% (95% CI, 32%–50%) lower risks of NSCL/P than GG genotypes in Asians, but these genotypes were not significant in Caucasians. For rs2013162, only AA was significant, that is, carried 0.65 (95% CI, 0.52–0.82) times lower odds than CC in Caucasians but not for Asians. For rs642961, AA and GA genotypes, respectively, carried 2.47 (95% CI, 1.41–4.35) and 1.40 (95% CI, 1.12–1.75) times higher odds in Asian, and 2.03 (95% CI, 1.52–2.71) and 1.58 (95% CI, 1.37–1.82) times higher odds in Caucasians compare with GG genotypes. For rs987525, AA and CA genotypes carried 2.27 (95% CI, 1.43–3.60) and 1.34 (95% CI, 1.02–1.77) times higher odds in Asian, and 5.25 (95% CI, 3.98–6.91) and 2.13 (95% CI–1.82, 2.49) times higher odds in Caucasians, and 1.42 (95% CI, 1.10–1.82) and 1.28 (95% CI, 1.09–1.50) times higher odds in mixed ethnicities compared with CC genotypes. These variant effects remained significant based on applying Bonferroni corrected‐thresholds, except in the mixed ethnicity. Conclusion We show robust variant effects in NSCL/P. Considering them with other genes and risk factors might be useful to improve prediction of NSCL/P occurrence. Birth Defects Research (Part A) 106:773–788, 2016. © 2016 The Authors Birth Defects Research Part A: Clinical and Molecular Teratology Published by Wiley Periodicals, Inc.


Introduction
Nonsyndromic cleft lip with or without cleft palate (NSCL/ P) is one of the commonest birth malformations and varies by geographic origin and ethnicity (Bender, 2000). The incidence is approximately 1 in 500 live births in Asia and America, 1 in 1000 in Europe, and 1 in 2500 in Africa. Cleft lip without cleft palate is more frequent in males than females, with a ratio of 2:1, and frequently occurs on the left rather than right side (Mossey et al., 2009).
Patients with NSCL/P suffer from functional and cosmetic problems, and NSCL/P can affect speech and communication, which results in delayed development, low self-confidence, and poor quality of life for both patient and family members. Although reconstruction with maxillo-facial surgery is possible, severe cases sometimes require multiple corrective surgeries (e.g., cheiloplasy, palatoplasty, rhinoplasty, orthognathicsurgery), dental work, and speech rehabilitation. Thus, this condition places a significant burden on family and medical services (Strauss, 1999;Murray, 2002).
There have been attempts to identify disease genes associated with NSCL/P as identification of disease genes may shed light on the etiology of the condition and facilitate efforts at prevention of disease. Mutations of interferon regulatory factor 6 (IRF6), located on chromosome 1q32.2 (OMIM # 607199), is associated with the autosomal-dominant Van der Woude syndrome, the most common Mendelian syndrome that has the cardinal signs of cleft lip with or without cleft palate (CL/P) and/or cleft palate only (CPO) with dental anomalies and pitted lips (Kondo et al., 2002;Rizos and Spyropoulos, 2004). Polymorphisms in IRF6 were later identified by Zucchero et al. (2004) and replicated by individual studies (Zucchero et al., 2004;Rahimov et al., 2008;Birnbaum et al., 2009b) and genome-wide association studies (Birnbaum et al., 2009b;Grant et al., 2009;Beaty et al., 2010;Ludwig et al., 2012), as associated with NSCL/P. Animal studies also indicate that IRF6 is involved in proliferationdifferentiation of the keratinocyte (Richardson et al., 2009), and hyper-proliferation of the epidermis may result in failure of terminal differentiation and multiple epithelial adhesion leading to CL/P (Ingraham et al., 2006;Richardson et al., 2006).
A genome-wide association study of German subjects of Central European origin identified 146 polymorphisms in the 8q24 region, of which rs987525 (OMIM # 612858) was the most significantly associated with NSCL/P and was the major susceptibility locus in the general population (Birnbaum et al., 2009b). Grant et al. (2009) also conducted a case-control genome-wide association studies (GWAS) in the United States and confirmed the associations of this region. In addition, Beaty et al. (2010) had replicated the GWAS in Asian subjects and had similar results, that is rs987525 was the most significant polymorphism. However, the biological mechanism of this SNP remains unclear.
Several independent studies replicated the effects of the IRF6 gene (at rs2235371 and rs642961) and the 8q24 region (at rs987525) on NSCL/P, and their effects were later pooled by a meta-analysis (Wang et al., 2012). However, this meta-analysis pooled only allele effects rather than genotypic effects, which assumes an additive mode of action. In addition, further studies have been published since then. We, therefore, performed an updated metaanalysis on the same polymorphisms with additional IRF6 polymorphisms at rs2013162 aiming to assess both variant effects and possible genetic models.

SEARCH STRATEGY
This review was conducted following the HuGE Review guideline Gwinn et al., 2014). Studies were identified by two reviewers (K.W. and S.R.) using Medline (by means of PubMed) and Scopus (Sci Verse Scopus; Elsevier B.V.) databases up to February 15, 2016. Search terms used were as follows: (irf6 or irf-6 or "interferon regulatory factor 6" or "interferon regulation factor-6" or 8q24) and (polymorphism or gene or allele or SNP) and (NSCL/P or "cleft lip" or CL or "cleft palate" or CP or CPO or "nonsyndromic cleft"). Details of the search strategies are described in Supplementary Table 1, which is available online.

INCLUSION AND EXCLUSION CRITERIA
Two reviewers (K.W. and S.R.) independently selected studies by looking through all titles and abstracts of the identified studies. Any human population-based (not familial-based) association study published in English was included if it met the following criteria: the outcome of interest was NSCL/P, the studied polymorphisms were in IRF6 or 8q24, and sufficient data were reported, that is, allele/genotype frequency between case and control; or odds ratio (OR) of allele/genotype effects on outcome along with its 95% confidence interval (CI). The studies with insufficient information were excluded if authors did not provide data after two contacts. Where there were multiple publications with the same or overlapping subjects (which could be detected by incorporating information on authorship lists, setting, study period, genotyping method, and subjects) the most complete and/or recent results were used. The reference lists of the selected articles were also reviewed to identify additional relevant publications.

DATA EXTRACTION
Summary data for IRF6 and 8q24 were extracted independently by two reviewers (K.W. and S.R.) using a standardized data extraction form. Data on co-variables, for example, mean age, percentage of males, ethnicity, and family history were also extracted. The main data (generic data and risk of bias assessment) were computerized separately by reviewer. Then, the two data were compared and validated to identify disagreement between. Meetings between the two reviewers (K.W. and S.R.) and the third party (A.T.) were then setup to discuss and solve disagreement.

RISK OF BIAS ASSESSMENT
The quality of studies was independently assessed by 2 reviewers (K.W. and S.R.) using a risk of bias assessment for genetic association studies, described in detail previously (Thakkinstian et al., 2005b(Thakkinstian et al., , 2011. Briefly, the assessment considered five domains: selection bias, information bias, population stratification, selective reports, and assessment of Hardy-Weinberg equilibrium (HWE). Each item was classified as yes, no, or unclear, corresponding to low risk of bias, high risk of bias, or unclear/insufficient information. Disagreement was solved by discussion with senior author (A.T.).

OUTCOME AND STUDIED GENE OF INTERESTS
Our review focused on NSCL/P as the outcomes of interest, which diagnosed according to original studied. The studied genes were three polymorphisms within the IRF6 gene (rs2235371 G>A, rs642961 G>A, and rs2013162 C>A), and one located within 8q24 (rs987525 C>A).

STATISTICAL ANALYSIS
Data analysis was performed following a standard method of meta-analysis for genetic association studies (Thakkinstian et al., 2005a(Thakkinstian et al., , 2012. HWE was checked in control groups for case-control studies using the Fisher's exact test. Data were then pooled based on studies that complied with HWE and if there were at least three studies for each studied polymorphism. A prevalence of minor allele frequency (MAF) for each polymorphism was estimated for each study, and then pooled across studies separately by ethnicity (Thakkinstian et al., 2005b).
Variant effects on NSCL/P were assessed using a perallele and per-genotype approach. For the per-allele approach, an OR of a versus A, if a and A were, respectively, minor and major alleles, along with its variance were estimated. Heterogeneity was then assessed using Cochrane's Q test and I 2 statistic. The ORs were then pooled using a random-effect model by Der Simonian & Laird method if heterogeneity was present (p value < 0.10 or I 2 25%); otherwise a fixed-effect model with inversevariance method was applied.
For the per-genotype approach, OR 1 (aa vs. AA) and OR 2 (Aa vs. AA) were estimated for each polymorphism across included studies. Heterogeneity was then assessed using the same methods described above. If any of the two ORs was heterogeneous, data were then pooled using a mixed-effect hierarchical logit model (Thakkinstian et al., 2012). In this model, the genotypes were considered as fixed effects, whereas the study was considered as a random effect. A likelihood ratio test was used to assess whether an overall variant effect was significant. The pooled ORs along with 95% CI were then estimated by exponential coefficients of the mixed-effect logit model. If the overall variant effect was present, the mode of genetic effect was captured by the parameter lambda (k 5 logOR 2 /logOR 1 ), which was then estimated using the model-free Bayesian approach (Minelli et al., 2005). The estimated k (range: 0 to 1) would suggests a recessive, dominant, and additive effect if k closes to 0, 1, and 0.5, respectively.
A sensitivity analysis was performed by including studies not observing HWE in the main pooling. Publication bias was assessed using the Egger's test and funnel plot (Egger et al., 1997). In addition, contour enhanced-funnel plots were generated if there was any evidence of asymmetry suggested by funnel plot or Egger's test (Peters et al., 2008).
Analyses were performed using STATA version 14 (Sta-taCorp, 2013;van Rooij et al., 2004) and WinBugs 1.4.2 (Spiegelhalter et al., 2007) with normal vague prior distributions for estimation of lambda and ORs. The model was run with 1000 iteration burn in, following by 10,000 iterations for estimation of parameters. A p-value less than 0.05 was considered statistically significant, except for tests of heterogeneity where a level of 0.10 was used. In addition, Bonferroni correction was applied to adjust for 27 multiple tests in total (see Supplementary Table 2), which resulted in a Bonferroni corrected threshold of 1.85 3 10 -3 .

IDENTIFYING STUDIES
A total of 307 and 1155 studies were, respectively, located from Medline by means of PubMed and Scopus (Fig. 1). After 282 duplicates were removed, 1180 studies were screened on titles or abstracts, and 34 studies were eligible. Reasons for ineligibility are provided in Figure 1, but were mainly nonrelevant genes or populations, and nongenetic association studies. For studies meeting inclusion criteria, the study by Letra et al. (2012) had insufficient data but the author provided additional data after our request. Studied polymorphisms and number of studies for each polymorphism are described in Figure 1. Among them, three polymorphisms (i.e., rs2235371, rs2013162, rs642961) for IRF6, and one polymorphism (i.e., rs987525) for 8q24 had sufficient numbers of studies; therefore, further pooling was focused on only these polymorphisms from 31 studies.

RISK OF BIAS ASSESSMENT
The results of bias assessment are presented in Supplementary Table 3. All studies had low risks of bias related to ascertainment of NSCL/P and population stratification. Twenty-nine of 31 studies assessed HWE. Most studies (97%) had low risk of selective reporting of outcomes. Twenty-seven (87%) and 24 (77%) studies described control ascertainment and quality control for genotyping, respectively.
Publication bias was assessed for both pooled ORs in Asian and Caucasian by funnel plot and Egger's test. There was no evidence of asymmetry of the funnel, see Supplementary Figure 1A-D. These results agreed with the Egger's tests, which indicated no evidence of asymmetry of funnels for both OR 1 (coefficient 5 1.448; p 5 0.408) and OR 2 (coefficient 5 -1.018; p 5 0.358) in Asian, and OR 1 (coefficient 5 0.480; p 5 0.261) and OR 2 (coefficient 5 -2.763; p 5 0.452) in Caucasian.
Publication bias was assessed for both pooled ORs in Caucasian and Asian by funnel plot and Egger's test. Funnel plots suggested little asymmetry (Supplementary Fig. 2A-D), but the Egger's test did not confirm this. The coefficients of asymmetry was -2.07 (SE 5 5.25; p 5 0.731) for OR 1 and -0.25 (SE 5 1.64; p 5 0.894) for OR 2 in Caucasian, and -8.58 (SE 5 5.45; p 5 0.360) for OR 1 and -7.41 (SE 5 1.88; p 5 0.158) for OR 2 in Asian.

Discussion
We conducted a systematic review and meta-analysis to assess genetic effects of polymorphisms in IRF6 and 8q24 region on NSCL/P. A total of 9 to 15 studies were included in the pooling for three polymorphisms (i.e., rs2235371, rs2013162, rs642961) for IRF6 and for rs987525 in 8q24. The MAF for these corresponding polymorphisms were similar across ethnicities except for rs2235371 where the MAF was common in Asians but rare in Caucasians, whereas MAF for rs987525 was common in Caucasians but less common in Asians. Because the MAF of rs2235371 was different between ethnicities, the genotype effect for this SNP was different. All polymorphisms were significantly associated with NSCL/P in some ethnicities. A significant protective effect of rs2235371 in Asians and effect of rs2013162 in Caucasians were suggested, that is, carrying one copy of the A allele for both polymorphisms lowered the odds of NSCL/P compared with G and C alleles, respectively. In addition, our findings indicated risk effects of minor alleles for rs642961 and rs987525. Our findings are consistent with previous human GWAS studies (Birnbaum et al., 2009b;Grant et al., 2009;Beaty et al., 2010;Ludwig et al., 2012). In addition, our findings are also consistent with the previous meta-analysis of genetic association studies by Wang et al. (2012), that included 10, 6, and 8 studies in the pooling of rs2235371 and rs642961 for IRF6, and rs987525 for 8q24, respectively. We confirmed the previous findings by expanding to a larger number of included genetic association studies (i.e., 15, 11, and 12 studies for these corresponding polymorphisms) and applying Bonferroni correction for multiple tests. Although each polymorphism was considered separately in meta-analysis for genetic association studies, this type of study had been used to replicate and confirm GWAS. Not only allele effects were estimated, but also pooled MAF and genotype effects by ethnicity additionally provided from previous evidences.
Formation of the lip and palate in human embryos starts at the fourth week and is complete by the twelfth week after fertilization. This process occurs by migration of the ectoderm, and formation of frontal, medial and lateral nasal processes, which then fuse together to form the normal lips and nose (Shkoukani et al., 2013). A normal palate also occurs after complete fusion of the primary and secondary palates; failure of fusion causes separation of the lip/s and/or palate, that is, cleft lip and cleft palate, respectively. Although the function of the IRF6 gene on the cleft lip formation is unclear, evidence has supported an effect of this gene on other birth malformations including Van der Woude syndrome and popliteal pterygium syndrome (Kondo et al., 2002;Kayano et al., 2003;Shotelersuk et al., 2003;Gatta et al., 2004;Ghassib e et al., 2004;Item et al., 2004). High IRF6 mRNA levels are also seen at the medial edge of the fusing palate (Kondo et al., 2002).
In addition, evidence from animal studies also supports the role of IRF6 on NSCL/P. Animal models have shown abnormalities of epithelial differentiation in mice with a mutation of IRF6 (Iwata et al., 2013). This supports a role for IRF6 in the formation and maintenance of the periderm, and spatiotemporal regulation for appropriate palatal adhesion (Richardson et al., 2009). Another study on mice with mutation of IRF6 showed abnormal skin, limb, and craniofacial development, which suggests a role for IRF6 in keratinocyte proliferation and differentiation of face and limb (Ingraham et al., 2006). The role of the 8q24 region on NSCL/P is still unclear but it might be involved in maintenance of an undifferentiated state in human acute myeloblastic leukemia cells (HL60) (Hirano et al., 2008), control of Myc expression (Uslu et al., 2014), or cranial facial enhancers (Attanasio et al., 2013).
The etiology of NSCL/P is still unclear but evidence suggests that both genetic and nongenetic factors play a role (Dixon et al., 2011;Shkoukani et al., 2013). Environmental factors may also play a role in NSCL/P. These factors may be independent or they may be modified by genes such as 5, 10methylenetetrahydrofolate reductase and folate . Evidence from cohort studies suggested that oral folate supplements before or during the periconception period might be able to reduce risk of NSCL/P (van Rooij et al., 2004;Kelly et al., 2012;Li et al., 2012); this is being confirmed by a randomized controlled trial (Wehby et al., 2012). Other factors associated with NSCL/P include alcohol consumption (Lieff et al., 1999;DeRoo et al., 2008), gestational diabetes (Correa et al., 2008), maternal age>40 years (Herkrath et al., 2012), and zinc deficiency (Warkany and Petering, 1972).
Our results also suggested that the IRF6 ((i.e., rs2235371, rs2013162, rs642961) and 8q24 (rs987525) SNPs contributed to NSCL/P, in which the rs987525 had highest contribution. This might be comparable to the effect of folate with an estimated PAR approximately 18% (van Rooij et al., 2004;Kelly et al., 2012). For public health prevention, screening these genes and other genes identified by previous meta-analyses (e.g., VAX1, Figueiredo et al., 2014; or bone morphogenetic protein 4, Hu et al., 2015), and combining them with environmental factors (e.g., folate, alcohol, etc.) may lead to a useful risk score for prediction of NSCL/P occurrence (Seddon et al., 2009).
Our study has some strengths. We attempted to pool variant effects on NSCL/P separately by ethnicity if data were available. Prevalence of MAFs for each study polymorphism was pooled across studies and ethnicity. Mode of variant effect was estimated using a mixed-effect logistic regression (Thakkinstian et al., 2005a) to capture lambda (Minelli et al., 2005). However, our study also has some limitations. First, given that we worked on summary data, we could not control for confounding effects, although the major source of confounding for genetic studies is population stratification and we summarized results by ethnic group. Second, we could only identify association between genes and NSCL/P, and we do not know whether the genes directly affect NSCL/P or whether their effects might be exerted through other environmental factors, such as folate. Exploring this would require individual patient data where data for folate and other covariables are also available. Third, we could not identify which polymorphism is the disease gene among those three polymorphisms in IRF6 gene because evidences showed linkage disequilibrium between rs2235371 and rs2013162 in Mexican (Ibarra-Arce et al., 2015) and Chinese (Mijiti et al., 2015) population (r 2 5 0.39 and 0.25, respectively) but not for rs642961 (Pan et al., 2010) (r 2 5 0.07). To determine this required raw data for performing a haplotype analysis. Some variant effects were heterogeneous (i.e., varied) across studies; we could explore for a cause of heterogeneity by fitting percent male in a meta-regression, but this was not detected (data were not shown). There might be other characteristics that cause heterogeneity, but those data were not available. Finally, there might be gene to gene interactions, but summary data do not allow us to explore this.

CONCLUSIONS
All polymorphisms, rs2235371, rs2013162, and rs642961 at IRF6 and rs987525 in 8q24 were significantly associated with NSCL/P. However, after applying Bonferroni corrected-thresholds, variant effects remained significant, except in mixed ethnicity. We show robust variant effects in NSCL/P. Including them with other risk factors and other genes in a risk prediction model of NSCL/P might be useful to improve prediction of NSCL/P occurrence.