Genetic Associations of Interleukin-related Genes with Graves’ Ophthalmopathy: a Systematic Review and Meta-analysis

Graves’ ophthalmopathy (GO) is the commonest extra-thyroidal manifestation of Graves’ disease (GD). Associations between interleukin-related (IL) gene polymorphisms and GO have been reported in different populations. We aim to confirm such associations by conducting a meta-analysis. Totally 382 publications were retrieved in MEDLINE and EMBASE up to 25/2/2015. After removing the duplicates and assessing the studies, we retrieved 16 studies that met the selection criteria for meta-analysis, involving 12 polymorphisms in 8 IL-related genes, and 1650 GO cases and 2909 GD controls. The summary odds ratio (OR) and 95% confidence intervals (CI) were estimated. We found one polymorphism in IL1A (rs1800587, c.-889C>T) showing a suggestive association with GO in the meta-analysis (allelic model [T vs. C]: OR = 1.62, 95% CI: 1.00–2.62, P = 0.050, I2 = 53.7%; recessive model [TT vs. TC + CC]: OR = 2.39, 95% CI: 1.07–5.37, P = 0.039, I2 = 23.6%; heterozygous model [TC vs. CC]: OR = 1.52, 95% CI: 1.04–2.22, P = 0.034, I2 = 37.0%). No association with GO was detected for the other 7 genes (IL1B, IL1RA, IL4, IL6, IL12B, IL13 and IL23R). Our results thus indicate that IL1A is likely to be a genetic biomarker for GO. Further studies with larger sample sizes are warranted to confirm the associations of IL1A and other IL-related genes with GO.

To explain the heterogeneity, we performed subgroup analysis by ethnicity. Due to the limited number of studies, we only tested the associations of 2 SNPs (IL1B rs16944 and IL4 rs2070874) in Chinese. However, these 2 SNPs did not show significant association with GO (P > 0.05), with low to moderate heterogeneities (Supplementary Table 2). Assessment of potential biases and sensitivity analysis. Lacka et al. compared a subgroup of patients with GD associated with GO from the onset and a subgroup contained patients in whom GO developed from 6 months to 7 years from the onset of GD 23 . To avoid selection bias, we conducted sensitivity analysis by excluding this study and keeping only patients with GD without GO as controls from the meta-analysis of IL1B rs16944 and rs1143634. The associations remained insignificant (Supplementary  Table 3). In the quality assessment of studies using the Newcastle Ottawa Scale (NOS), all of the studies were assigned 7 or more stars, indicating low risk of introducing biases. Therefore, no study was excluded from the meta-analysis due to poor quality (Supplementary Table 4 more studies (i.e., IL1B rs1143634 and rs16944, IL4 rs2070874 and IL13 rs1800925), sensitivity analyses were performed by sequentially omitting one study at a time. The insignificant associations remained unchanged (P > 0.05; data not shown). There was no significant publication bias detected by the funnel plots (data not shown) and Egger's test (Table 2 and Supplementary Table 3).

Discussion
This study has, for the first time, summarized the associations of IL-related genes with GO. Among the 11 reported genes, we performed meta-analyses on 12 polymorphisms in 8 genes. Unexpectedly, we found only one SNP, IL1A rs1800587 (c.-889C> T), being marginally associated with GO. No significant association was detected for SNPs in the other 7 genes (IL1B, IL1RA, IL4, IL6, IL12B, IL13 and IL23R), among which IL1B rs1143634, IL1RA A2/non-A2 and IL23R rs7530511 showed no heterogeneity across the study populations. The IL1A SNP rs1800587 showed a suggestive association with no to moderate heterogeneities in different genetic models. This SNP was reported in 2 studies 21,22 . Although a significant association was reported only in one study 21 , the effect of the risk allele T pointed to the same direction in the both studies (OR = 2.16 21 and OR = 1.32 22 ). The heterogeneity could be due to the relatively small sample size in the study of Khalilzadeh et al. (about 50 cases and 50 controls), the ethnic differences in linkage disequilibrium structures, and the differences in the minor allele frequencies (MAF) of rs1800587 (T) between Iranian (about 43.0%) 21 and Chinese (about 10%) 22 . Notably, however, since the P values did not survive the Bonferroni correction for multiple testing, the genetic association of the IL1A SNP with GO has yet to be confirmed in further studies with larger sample sizes.
IL1A, a major member of the IL1 superfamily, is the prototype pro-inflammatory and a potent pleiotropic cytokine involved in acute or chronic inflammation 50 . Associations between IL1A and GO were demonstrated in biochemical, histological, immunological and genetics studies. There were significant differences in the serum IL1A levels between controls and GO patients, and for the latter, before and after corticosteroid, corticosteroid with orbital irradiation, or decompression 51 . IL1A immunoreactivity was detected in the orbital tissues, their fibroblast cultures and supernatants from 5 out of 6 GO patients, but absent in those derived from 5 normal individuals 52 . An in vitro study demonstrated the induction of intercellular adhesion molecule 1 (ICAM-1), endothelial leukocyte adhesion molecule 1 (ELAM-1) and vascular cell adhesion molecule 1 (VCAM-1), which promote T cell chemotaxis upon the exposure of endothelial cells generated from retrobulbar tissues to IL1A 53 . The proliferation of orbital fibroblasts from GO patients was stimulated by IL1A, which has no effect on normal orbital fibroblasts 54 . Transcription of prostaglandin endoperoxidase H synthase-2 (an inflammatory cyclooxygenase that produces prostaglandin E2 and contributes to orbital inflammation in GO 55 ) in orbital fibroblasts by leukoregulin (a product of activated T lymphocytes) was found to be mediated through an intermediate induction of IL1A 56 . In our meta-analysis, we found the IL1A SNP rs1800587 as a potential susceptibility genetic marker for GO, confirming the involvement of IL1A in the disease. In fact, the IL1A SNP rs1800587, located in the 5' untranslated region (c.-889C> T), had been associated with autoimmune diseases including ankylosing spondylitis 57 , systemic lupus erythematosus 58 , psoriatic arthritis 59 and Behcet's disease 60 . IL1A and its SNP rs1800587 could thus play a role in the pathogenesis of autoimmune diseases including GO.
Except for IL1A, SNPs in other reported genes did not show a significant association with GO in our meta-analysis. Three SNPs showed no association with GO in all of the tested populations (P > 0.05) with no heterogeneity, including IL1B rs1143634 (Caucasians, Iranian and Chinese), IL1RA A2/non-A2 (Caucasians) and IL23R rs7530511 (Caucasians and Japanese). They are not likely to be genetic markers for GO. Another 6 insignificant SNPs (IL1B rs16944, IL4 rs2070874, IL6 rs1800795, IL13 rs1800925 and c.-2044G> A, and IL23R rs10889677) also lacked significant association in any of the studies with mild to high heterogeneities. In contrast, IL12B rs3212227 28 and IL23R rs2201841 41 showed significant associations with GO in Iranian and Caucasians, respectively, but not in Japanese 34,42 . In our meta-analysis, no significant association was found for these 2 SNPs by using the random-effect model, with moderate to high heterogeneities. Of note, the IL23R SNP rs2201841 was also significantly associated with other autoimmune diseases, such as Crohn's Disease 61,62 and rheumatoid arthritis 62 . Therefore, further replication of these 2 SNPs, IL12B rs3212227 and IL23R rs2201841, in GO among different populations are warranted.
Two IL1B SNPs, rs1143634 and rs16944, showed no association with GO. The summary results of the IL1B gene in our meta-analysis were inconsistent with that in the study of Liu et al. 22 . In this study, we used patients with GD but without GO as controls, with a view to assess the effects of the gene SNPs on GO in a background of GD. In contrast, Liu's group compared GO patients with healthy subjects and detected a significant association 22 . Thus, our results cannot be compared directly with that of Liu et al. Further studies are needed to confirm if the IL1B SNPs are genuine markers differentiating GO patients from normal subjects.
IL1RA acts as a competitive inhibitor of IL1A and IL1B and blocks IL1-mediated cellular activities 63 , such as IL-1-induced glycosaminoglycan production by cultured human orbital fibroblasts 64 . IL1RA could also block the induction of prostaglandin endoperoxidase H synthase-2 by leukoregulin 56 . A study had shown that upon cytokines exposure, markedly lower level of IL1RA expression was found in cultured orbital fibroblasts of GO patients as compared to the normal orbital fibroblasts 65 . Although IL1RA A2/non-A2 was not significant in our meta-analysis, another IL1RA SNP (c.11100C> T) had shown positive association in Iranians 21 . Follow-up studies on the IL1RA polymorphisms are needed to confirm the role of IL1RA in GO in specific population such as Iranians.
IL4 is a potent Th2 cytokine which stimulates proliferation of IgE-and IgG-secreting B cells and the expression of HLA class II antigens via STAT6 66 against Th1 inflammatory response 30 . It has also been detected in orbital fat tissues of GO patients 12 . The promoter SNP rs2070874, which has transcriptional activity 30,67 , did not show a significant association with GO in our meta-analysis. However, significant Scientific RepoRts | 5:16672 | DOI: 10.1038/srep16672 associations of other IL4 SNPs, including c.-1098T> G and c.-33C> T, with GO have been reported in Iranians 31 . Therefore, these two IL4 polymorphisms should be tested in future studies. IL13 and IL4 have similar biological functions 68 . IL13 is an anti-inflammatory cytokine that regulates IgE synthesis 68,69 and the maturation of B cells 30 . However, there was no significant association detected in Chinese (rs1800925) 30 , Japanese (rs1800925 and c.-2044G> A) 37 and Caucasians (rs1800925 and c.-2044G> A) 36 . Consistently, we did not detect associations of these two IL13 SNPs with GO.
The IL6 SNP rs1800795 is associated with multiple autoimmune diseases, including systemic-onset juvenile chronic arthritis 70 , type I diabetes mellitus 71 , rheumatoid arthritis 72 and Sjogren's syndrome 73 . When compared with healthy controls, serum IL6 levels were significantly higher in GD and GO patients, especially in active GO patients 74 . However, only one SNP rs1800795 in IL6 was eligible for the meta-analysis, and it showed a lack of significant association.
This meta-analysis also reveals several limitations in the existing genetic studies of GO. First, the small number of published genetic studies on IL-related genes in GO limited the power of determining the associations, especially among different ethnic groups. Second, GO may not develop concurrently with GD. Classifying GO based on a cross-sectional assessment of observer-dependent signs and subject-dependent symptoms may therefore introduce bias. Third, as the pathogenesis of GO and GD is multifactorial, it would be more informative to test genetic, environmental (e.g. smoking), hormonal (e.g. fluctuation of thyroid function) and antigenic (thyroid related autoantibodies and use of radioactive iodine) factors, and their interactions in the study population. However, few genetic studies on IL-related genes and GO provided such information.
In conclusion, in this systematic review and meta-analysis of the association of IL-related genes with GO, we identified IL1A rs1800587 as the only SNP that is potentially associated with GO. Since the overall number of studies is small, further studies with larger sample sizes are needed to confirm IL1A rs1800587 as a genetic biomarker for GO, and also verify the roles of other IL-related genes in the disease.

Methods
Searching Strategy. We performed the literature search in the EMBASE and MEDLINE databases via the Ovid platform using structured search strategies. We identified citations recorded during the period starting from May 1, 1989 to February 25, 2015. Boolean logic and search terms with controlled vocabularies (i.e. Medical Subject Heading terms) were used: (Graves' ophthalmopathy OR thyroid associated ophthalmopathy) AND interleukins (Supplementary Table 5). Moreover, we manually scanned the reference lists of the identified articles, reviews, and meta-analyses to include all potentially relevant articles. No language filters were applied in the literature search.

Inclusion and Exclusion Criteria.
A study was included if it fulfilled the following criteria: (1) original case-control study on the associations of IL-related genes polymorphisms with GO; (2) cases were patients with GO defined by clinical signs of GO or NOSPECS Class 2 or 3 and above; (3) controls were patients with Graves' disease without GO (no clinical signs of GO or NOSPECS Class 0 or 1); (4) study subjects were unrelated individuals from clearly defined populations; (5) allele or genotype counts or frequencies in both case and control groups were provided (or existing data allow their calculation). Animal studies, case reports, reviews, abstracts, conference proceedings, editorials and studies with incomplete data were excluded. Literature Review and Data Extraction. Two investigators (W.K.H. and S.S.R.) screened and reviewed all studies independently. Disagreement was resolved by thorough discussion with a third investigator (L.J.C.) until consensus was reached. A customized data form was used to extract the data, which included the first author, year of publication, country of study, ethnicity, definition of cases and controls, sample size in case and control groups, gene and polymorphisms studied, allelic and genotypic counts, and result of the test for HWE in the control group. Two reviewers (W.K.H. and S.S.R.) extracted the data independently. Disagreement was resolved by consensus among the investigators. If the allele counts were not reported, we calculated them from the genotype data. If genotype counts were missing, we estimated the data using allele frequencies (if available) and sample sizes, assuming no deviation from HWE unless otherwise reported 75 . If there was no extractable genetic information in an eligible study, we communicated with the authors for the data. Allele counts of the eligible SNPs for meta-analysis were summarized in Supplementary Table 6.
Statistical Analysis. Meta-analysis for each gene polymorphism was performed if it was reported in 2 or more studies. The genetic association was assessed using different genetic models, including allelic (A vs. a), dominant (AA+Aa vs. aa), recessive (aa vs. AA+Aa) and codominant (homozygous: AA vs. aa; heterozygous: AA vs. Aa) models. The strength of association was evaluated using the summary odds ratios and 95% confidence intervals of each gene polymorphism. Heterogeneity was tested by the Q-statistic and the I 2 value 76 . The Q-statistic was considered significant when P < 0.10. The I 2 values indicated no (0-24.9%), low (25-49.9%), moderate (50-74.9%) or high (75-100%) inter-study heterogeneity 76,77 . If the P value for the Q-statistic was < 0.10 or the I 2 value ≥ 50%, a random-effect model was used, otherwise a fixed-effect model was adopted 78 . In the assessment of data quality, we first examined the HWE in the control groups. If HWE was not reported, we tested it using the control group data with the Chi-square test. Also, we adopted the Newcastle Ottawa Scale (accessed via http://www.ohri.ca/programs/clinical_epidemiology/oxford.asp) to evaluate the quality of case-control studies (Appendix 1). A star was given to each study if one requirement in the NOS from 3 dimensions (selection, comparability and exposure) was met. The maximum number of stars that can be assigned to one study was 9. A study obtaining ≤ 6 stars was considered as of high risk in introducing bias 79 . We then conducted a sensitivity analysis to confirm the associations by sequentially omitting each of the studies one at a time, studies deviated from HWE, or studies of suboptimal quality 80 . Furthermore, the Funnel plots and Egger's test were performed to assess potential biases (e.g. publication bias) 81,82 . The presence of bias was suggested when the P value of the Egger's test was < 0.05. All statistical analyses were performed using the R software for statistical computing (v3.0.0, http://cran.r-project.org/). Of note, since we tested the genetic association using 5 genetic models, the risk of type I error might be increased; therefore, we corrected the P values for association using the Bonferroni method. Thus, a P value of less than 0.010 (0.05/5) was considered statistically significant, where 5 is the number of genetic models being tested.