Abstract
Background: Linkage analyses and association studies suggested that inherited genetic variations play a role in the development of differentiated thyroid carcinoma (DTC).
Methods: We combined the results from a genome-wide association study (GWAS) performed by our group and from published studies on DTC. With a first approach, we evaluated whether a SNP published as associated with the risk of DTC could replicate in our GWAS (using FDR as adjustment for multiple comparisons). With the second approach, meta-analyses were performed between literature and GWAS when both sources suggested an association, increasing the statistical power of the analysis.
Results: rs1799814 (CYP1A1), rs1121980 (FTO), and 3 SNPs within 9q22 (rs965513, rs7048394, and rs894673) replicated the associations described in the literature. In addition, the meta-analyses between literature and GWAS revealed 10 more SNPs within 9q22, six within FTO, two within SOD1, and single variations within HUS1, WDR3, UGT2B7, ALOX12, TICAM1, ATG16L1, HDAC4, PIK3CA, SULF1, IL11RA, VEGFA, and 1p31.3, 2q35, 8p12, and 14q13.
Conclusion: This analysis confirmed several published risk loci that could be involved in DTC predisposition.
Impact: These findings provide evidence for the role of germline variants in DTC etiology and are consistent with a polygenic model of the disease. Cancer Epidemiol Biomarkers Prev; 25(4); 700–13. ©2016 AACR.
Introduction
Thyroid carcinoma is the most common malignancy of the endocrine system showing increasing incidences over the years (1–3), with age-standardized rates (ASR) of about 6/100,000 in the developed countries. Particularly elevated ASRs were observed in Lithuania (ASR = 15.5/100,000), Italy (ASR = 13.5/100,000), Austria (ASR = 12.4/100,000), and in the United States (ASR = 9.9/100,000; refs. 4, 5). Two main thyroid carcinoma histological types can be distinguished: the “medullary” and the “non-medullary” thyroid carcinoma, the former (MTC) originating from the para-follicular cells, the latter (NMTC) from the follicular cells. NMTC comprises the most frequent subtypes, papillary (PTC) and follicular (FTC) thyroid carcinomas (defined overall as “differentiated thyroid carcinomas”, DTC), accounting for 80% and 15% of the cases, respectively. Hürthle cells (or oxyphilic cells, 5%) and poorly differentiated carcinomas (1%–6%) are considered as not common entities (6, 7).
The great majority of DTCs behaves as a sporadic form, featured by somatic mutations within RET, RAS, BRAF, or NTRK1 genes and affecting the MAPK signaling pathway (8–11). However, approximately 5% of cases, mostly PTC, have a family history (12). Inherited genetic variations play an important role in both the familial and the sporadic forms, as supported by data from linkage analyses or case–control association studies (CCASs). In particular, CCASs, when carried out with an appropriate sample size, constitute the state-of-the-art for the identification of common genetic variants associated with complex traits (13). In the literature, a large number of gene-based studies have been published in which specific a priori hypotheses have been examined. Thus, SNPs within genes encoding for proteins involved in the DNA repair, cell-cycle control, thyroid physiology, or playing a role in other types of human cancer have been investigated (14). The major weakness of these studies is that unknown genes playing a relevant role in the etiology of the disease could have been missed. This limitation is solved by genome-wide association studies (GWASs) in which the whole genome is analyzed without formulating any a priori hypothesis. GWASs on DTC allowed discovering novel variants, including those near FOXE1, DIRC3, NKX2-1 (15–17) and, more recently, those near IMMP2L, RARRES1, SNAPC4/CARD9, ARSB, BATF, DHX35, SPATA13, GALNTL4, and FOXA2 (18–20). However, to ensure a high quality and to prevent false-positive findings, highly stringent criteria are applied in the GWASs with the disadvantage of excluding SNPs truly associated with the risk. In the present work, we investigated whether SNPs associated with the susceptibility to DTC in previous CCASs could replicate in an independent GWAS carried out by our research group. Moreover, we investigated whether SNPs showing a sub-threshold genome-wide statistical significance in our GWAS could improve their association following a meta-analysis with previously published data.
Materials and Methods
Ethics statement
Study participants were recruited according to the protocols approved by the institutional review boards in accordance with the Declaration of Helsinki. All subjects provided written informed consent to participate in the study and allowed the use of their biological samples.
Study participants of the GWAS
The group of cases comprised 701 histologically confirmed DTC patients from central and southern Italy, recruited at the Cisanello Hospital in Pisa, an important Italian referral center for thyroid diseases. The control group comprised 499 healthy individuals from the Meyer Hospital in Florence without known thyroid disease, of which 390 were blood donor volunteers and 109 were healthy individuals recruited during a routine health screening. Cases and controls were frequency matched by sex, age, body mass index (BMI), and smoking habits. The patient group consisted of 22.3% males and 77.7% females with a median age of 46; the control group consisted of 23.2% males and 76.8% females with a median age of 50. The median BMI was 24.5 in cases and 24.4 in controls. The proportion of smokers was 37% in cases and 40% in controls. All cases and controls were of Caucasian origin.
Genome-wide association study
Full details of the GWAS, including the genotyping process, quality control and statistical analysis were previously described (20). Briefly, samples were genotyped using Illumina HumanOmni1-Quad_v1-0_B 1M BeadChips and Illumina HumanOmniExpress-12v1_A 730K BeadChips. Genotype calling was performed using Illumina GenomeStudio 2010 (Illumina Inc.). After applying strict quality control criteria, the analysis was restricted to the subset of genotyped SNPs common to both Illumina arrays used. Hence, 572 042 SNPs were analyzed for association with DTC risk in 690 cases and 497 controls. The adequacy of the case–control matching and the possibility of differential genotyping of cases and controls were assessed using Q-Q plots of test statistics. The genomic control inflation factor λ was calculated using the standard method by the CRAN R package GAP (Genetic Analysis Package; https://cran.r-project.org/web/packages/gap/index.html; http://www.inside-r.org/packages/cran/gap/docs/gcontrol2).
The inflation factor λ was 1.0, excluding the possibility of hidden population substructure, relatedness among subjects or differential genotype calling. Statistical analysis was conducted using PLINK version 1.06 (21).
Search strategy and selection criteria
PubMed was searched from database inception until September 2013 to collect case–control studies investigating the association between SNPs and DTC. We used the keywords polymorph* AND (papillary OR follicular OR non-medullary OR “non medullary”) AND thyroid AND (cancer OR carcinoma) AND (susceptibility OR risk OR predisposition) to collect studies carried out on DTC or PTC. The major reasons for exclusion of the studies were (i) studies not in English language; (ii) studies without odds ratio (OR) and 95% confidence interval (95% CI); (iii) case–case studies; (iv) studies on benign thyroid disease. A total of 100 original articles and five meta-analyses met our criteria and were assessed. The list of citations is reported, for brevity, in the Supplementary data. The SNPs reported in these studies were recorded and searched in the present GWAS, allowing a direct comparison between the results published in the literature with the results from the GWAS. When an SNP was not found in the GWAS, the linkage disequilibrium (LD) block around the SNP was checked using the CEU data of the 1000 Genome Project (22), and the results of SNPs in high LD (r2 ≥ 0.8) were reported. All the collected data are reported in Supplementary Tables S1 and S2.
Statistical approaches
Two statistical approaches were used to reduce the number of false positives and to increase the power of the study. With the first approach, we performed a meta-analysis of published data when more than one study was carried out on a given SNP. Then, we evaluated whether SNPs previously associated with the risk of DTC (positive, at a nominal significance level of Pass < 0.05, either in a single study or in meta-analysis) were associated also with the risk of DTC in our GWAS. These SNPs were evaluated by calculating their allelic Pass in the GWAS. In order to adjust for multiple comparisons, the false-positive discovery rate correction (FDR; ref. 23) was applied to the list of Pass obtained in the GWAS and the associations with q < 0.05 were considered as statistically significant, i.e., considered as replicating the literature data.
With the second approach, results from SNPs positive in the literature (either in a single study or in the meta-analysis of the literature) were meta-analyzed with those of GWAS. Moreover, the meta-analyses were performed also when, for a given SNP, a suggestive evidence of association (Pass < 0.20, taken arbitrarily) was observed both in the literature and in the GWAS. Because the GWAS was performed on Caucasians, the meta-analysis first was carried out in Caucasians. When literature data were not available for Caucasians, the GWAS was meta-analyzed using literature data for the available population(s).
The pooled ORs were calculated for allelic model (a vs. A) and additive model (Aa vs. AA and aa vs. AA). In case only dominant or recessive model was reported in the literature, the same model was applied for the GWAS data.
The statistics are based on the absolute counts of variant and common alleles/genotypes among cases and controls. The χ2 based Q-test was used to assess heterogeneity across studies (Phet < 0.05) and I2 statistics was calculated to quantify the proportion of the total variation across studies due to heterogeneity. In case of no significant heterogeneity, OR and 95% CI were assessed using the fixed-effect model (the Mantel–Haenszel method); otherwise the random-effects model (DerSimonian–Laird method) was used. Meta-analyses were performed by MIX 1.7 freeware software. Also in this case, adjustment for multiple comparisons was performed by applying the FDR correction and q < 0.05 were considered as significant. A SNP associated with the risk of DTC in the literature was considered replicated when found with a q < 0.05 also in the GWAS. Moreover, an SNP was considered positively associated with the risk of DTC when found with a q < 0.05 in the meta-analysis.
Results
One hundred published articles, reporting results for 316 SNPs belonging to 127 genes, met the selection criteria (see the reference list in the Supplementary data). Data collected included the reference of the literature, the gene name, the dbSNP identification number, the number of cases and controls investigated, and the OR with its 95% CI, of the allelic and additive models. In the first type of evaluation, the corresponding ORs and 95% CIs were also calculated for these SNPs based on the GWAS. In case only dominant or recessive model was reported in the literature, the same model was applied for the GWAS data. The results are reported side-by-side in Supplementary Tables S1 and S2, respectively, to allow a direct comparison. Among the 316 SNPs, 91 were associated with the risk of DTC in a statistically significant way according to the literature (Pass < 0.05). The meta-analysis of the literature data alone was performed on 46 SNPs and 13 were statistically significant at the 0.05 level (Supplementary Table S3). Fifteen of the 91 SNPs associated in any study were replicated in the GWAS at the same significance level, and the side-by-side comparison is shown in Table 1. However, only five SNPs, including CYP1A1 rs1799814, FTO rs1121980, and the GWAS identified SNPs on 9q22 (rs965513, rs7048394, and rs894673), were statistically significant after the application of FDR correction. Only one SNP (rs965513) showed to be associated in a statistically significant way in the meta-analysis of literature data and in the present GWAS (Table 1). In addition to these analyses, we adopted another approach. In order to ascertain whether an increase of statistical power could allow reaching a statistical significance, we selected SNPs showing a Pass < 0.2 (arbitrarily chosen, in any inheritance model) both in the literature and in the GWAS and we performed a meta-analysis. Moreover, we added the results of the present GWAS to meta-analyses from literature data, when these latter showed SNPs significantly associated with the risk of DTC.
List of SNPs associated with the risk of DTC in previous studies published in the literature
Results from the meta-analyses are shown in Tables 2 and 3. Below, we summarize the results according to the biological processes.
Meta-analyses of published data on Caucasians with data from present GWAS
Meta-analyses between the non-Caucasian population(s) from the literature and present GWAS
SNPs within DNA repair genes
A total of 64 SNPs located within 27 genes involved in DNA repair pathways were investigated so far in the context of DTC. Of them, 10 were associated with the disease in the literature with at least one genetic model but none was replicated in the present GWAS under the allelic model (Table 1). However, a statistical significance was observed for rs25487 within XRCC1 (OR = 0.76, 95% CI, 0.59–0.99 for heterozygotes in the additive model). When this result was combined in meta-analysis with seven previous studies carried out on Caucasians, an OR of 0.92 (95% CI, 0.85–0.99) was found in the allelic model; however, this result was not significant after FDR correction (q = 0.06). Moreover, the meta-analyses revealed an increased risk for rs2708906, at 5′ region near HUS1 (OR = 1.34; 95% CI, 1.08–1.64; q = 0.04 for heterozygotes and OR = 1.52; 95% CI, 1.16–2.00; q = 0.01 for homozygotes; Table 2).
SNPs within cell-cycle regulation and apoptosis genes
Thirty-three common SNPs in 15 genes involved in the cell-cycle regulation or in apoptosis were collected. A total of nine significance associations were published so far, but only rs4658973 (WDR3) was replicated in the GWAS (allelic model: OR = 0.83; 95% CI, 0.70–1.00). Meta-analysis on Caucasians again suggested the role of this variant in DTC etiology with OR = 0.71 (95% CI, 0.61–0.82, allelic model), remaining significant after FDR correction (q = 1.8 × 10−6; Tables 1 and 2). Moreover, when the GWAS association on Caucasians was combined with a mixed population from a previous study, a statistical significance was found for rs2279744 (MDM2, OR = 1.40, 95% CI, 1.12–1.74; only recessive model was available for meta-analysis; Table 3).
SNPs within genes encoding for xenobiotic metabolism enzymes
Through PubMed search, 67 SNPs within 19 genes encoding for xenobiotic metabolism enzymes (XME) were collected. Overall, 19 positive associations were reported in the literature. Interestingly, rs1799814 (CYP1A1) showed a strong association in GWAS and it remained statistically significant after FDR correction (Table 1). In the meta-analysis, a high risk was found associated with the rare allele (OR = 1.86, 95% CI, 1.50–2.30, q = 4.4 × 10−8). Besides rs1799814, meta-analyses on Caucasians, after multiple testing correction, revealed a possible role also for rs1041740 (SOD1, q = 5.5 × 10−3; allelic model), rs12626475 (3′ region near SOD1, q = 0.02; allelic model), and rs3924194 (UGT2B7, q = 0.04; for heterozygotes; Table 2).
SNPs within genes involved in thyroid function
Seven genes playing a key role in thyroid function were assessed in DTC studies by genotyping 21 SNPs. Only five SNPs were reported as significantly associated with DTC, and none of them was significant in the GWAS. Thus, the meta-analyses did not confirm the role of these variants in DTC etiology (Tables 1 and 2).
SNPs within MAPK pathway genes
Of 17 SNPs within 8 genes of the MAPK pathway, four predisposing variants were reported. GWAS replicated the significant association found for rs12628 (HRAS, OR = 1.23, 95% CI, 1.02–1.48 in the allelic model), but a high heterogeneity was found between the study population previously analyzed and the present study (Phet < 0.0001). Thus, no significant evidence of association was identified in the meta-analysis using the random-effect model (Table 2).
SNPs within immune response and inflammation genes
Fifteen genes and 33 SNPs involved in immunity or in inflammation pathways were analyzed to identify susceptibility variants for DTC and eight significantly associated SNPs were published (Table 1). The present GWAS replicated the possible role of rs1126667 (ALOX12, OR = 1.34, 95% CI, 1.02–1.75; heterozygotes) and rs2292151 (TICAM1, OR = 1.69, 95% CI, 0.99–2.91; homozygotes; Table 2), as well as for rs1061758 (IL11RA, OR = 1.29, 95% CI, 1.02–1.65; allelic model; Table 3). The involvement of these SNPs in increasing the risk of DTC was further suggested by the meta-analyses. In particular, an OR = 1.74 (95% CI, 1.28–2.37, q = 8.9 × 10−4; allelic model) for rs1126667, and OR = 1.24 (95% CI, 1.06–1.45, q = 0.01; allelic model) for rs2292151 was observed in the meta-analysis with Caucasians studies, and an OR = 1.39 (95% CI, 1.14–1.70, q= 0.01; allelic model) for 1061758 was observed in the meta-analysis with an Asian study (Tables 2 and 3).
SNPs within other cancer genes
Fifty-five SNPs in 28 other genes related to cancer were investigated in relation to DTC risk and 37 SNPs were associated according to the literature. Of them, SNPs within ATG16L1 and FTO showed a strong association in GWAS under the allelic model and rs1121980 (within FTO) remained associated also after FDR correction (Table 1). According to the present meta-analysis of the published results and our GWAS data on Caucasians, SNPs rs2241880 (ATG16L1, OR = 0.81, 95% CI, 0.70–0.93, q = 7.6 × 10−3; allelic model), rs11642841 (FTO, OR = 0.76, 95% CI, 0.67–0.87, q = 1.2 × 10−4; allelic model), rs1121980 (FTO, OR = 0.75, 95% CI, 0.66–0.86, q = 5.7 × 10−5; allelic model), rs8050136 (FTO, OR = 0.76, 95% CI, 0.67–0.86, q = 4.8 × 10−5; allelic model), rs9939609 (FTO, OR = 0.77, 95% CI, 0.67–0.88, q = 3.9 × 10−4; allelic model), rs7202116 (FTO, OR = 0.76, 95% CI, 0.66–0.87, q = 2.5 × 10−4; allelic model), rs7584828 (HDAC4, OR = 0.68, 95% CI, 0.54–0.84, q = 4.2 × 10−3; heterozygotes), rs2132572 (5′ region near IGFBP3, OR = 0.77, 95% CI, 0.61–0.96, only dominant model was available for meta-analysis), and rs17849071 (PIK3CA, OR = 0.64, 95% CI, 0.46–0.90, q = 0.04; heterozygotes) were associated with a reduced risk of DTC, whereas SNPs rs17817288 (FTO, OR = 1.32, 95% CI, 1.15–1.51, q = 1.6 × 10−4; allelic model) and rs6472462 (5′ region near SULF1, OR = 1.17, 95% CI, 1.03–1.33, q = 0.03; allelic model) were associated with increased risks (Table 2). When the meta-analyses were extended to other available populations, four more SNPs showed an evidence of association, although not significantly after FDR correction: rs2229765 (IGF1R, OR = 0.77, 95% CI, 0.61–0.98, q = 0.09; heterozygotes), rs2230396 (ITGB1, OR = 0.75, 95% CI, 0.58–0.98, q = 0.09; heterozygotes), rs17524488 (5′ region near OPN, OR = 0.81, 95% CI, 0.67–0.99, q = 0.09; heterozygotes), and rs699947 (5′ region near VEGFA, OR = 1.22, 95% CI, 1.05–1.41, q = 0.05; allelic model; Table 3).
SNPs previously studied in relation to DTC risks from genome-wide association studies or studies focused on specific intergenic regions
Genetic variants on 1p31.3, 2q35, 8p12, 9q22, and 14q13.3 were associated with DTC risk by using genome-wide approaches. Three LD blocks (defined by rs965513, rs7048394, and rs894673) on chromosome 9q22 near FOXE1 were associated with DTC risk so far. These SNPs also showed a strong association in the present GWAS, where the allelic Pass remained statistically significant also after FDR correction (Table 1). Moreover, these associations were strengthened in the meta-analyses on Caucasians in the allelic model, with OR = 1.85 (95% CI, 1.76–1.95, q < 10−20) for rs965513, OR = 1.51 (95% CI, 1.31–1.73, q = 2.3 × 10−8) for rs7048394 and OR = 1.51 (95% CI, 1.33–1.71, q = 8.3 × 10−10) for rs894673. Moreover, the present meta-analysis points rs334725 (1p13.3, OR = 1.32, 95% CI, 1.10–1.59, q = 5.1 × 10−3), rs966423 (2q35, OR = 1.27, 95% CI, 1.19–1.35, q = 1.3 × 10−12), rs2439302 (8p12, OR = 1.30, 95% CI, 1.23–1.39, q = 1.2 × 10−15), and rs944289 (14q13, OR = 1.25, 95% CI, 1.17–1.33, q = 0.02) as associated with the risk of DTC (Table 2).
Discussion
Current scientists' knowledge on DTC genetic risk factors is based on a series of association studies on genes involved in different cellular mechanisms that could lead to malignant transformation of thyroid cells. Typically these studies were performed according to candidate-gene approaches, and rarely the findings were replicated using similar samples in terms of ethnicity and thyroid carcinoma histological type. Furthermore, to date, only few GWASs were performed, and a small number of genomic loci were associated with the risk of the disease by using this approach.
In order to gain further insights into the role of SNPs previously associated with DTC, in the present work we carefully analyzed the results of our GWAS and we performed meta-analyses with the previous studies. The associations between DTC and well-established GWAS-identified SNPs, including rs965513, rs7048394, and rs894673 near FOXE1 (9q22), were replicated using our GWAS data. Furthermore, rs944289 near NKX2-1 (14q13.3), rs966423 within DIRC3 (2q35), rs334725 within NFIA (1p31.3), and rs2439302 within NRG1 (8p12) showed an evidence of association in the meta-analysis of the GWAS results and previous published data. The role of these loci in DTC etiology was already discussed in previous works and will be not discussed here.
Although in the present work most of the SNPs assayed in previously published hypothesis-driven studies were not associated with the risk of DTC, it is noteworthy to observe that several of them actually did associate. In particular, rs1799814 within CYP1A1 and rs1121980 within FTO were replicated on our GWAS data after the application of multiple testing corrections. The meta-analysis–based approach provided an evidence of association of several additional variants, including SNPs in the DNA repair gene HUS1 (rs2708906), cell-cycle regulation gene WDR3 (rs4658973), xenobiotic metabolism genes SOD1 (rs1041740, rs12626475) and UGT2B7 (rs3924194), the immune response and inflammation genes ALOX12 (rs1126667), TICAM1 (rs2292151), and IL11RA (rs1061758), as well as other cancer genes ATG16L1 (rs2241880), FTO (rs17817288, rs11642841, rs9939609), HDAC4 (rs7584828), IGFBP3 (rs2132572), PIK3CA (rs17849071), SULF1 (rs6472462), IGF1R (rs2229765), OPN (rs17524488), and VEGFA (rs699947). All these SNPs were previously investigated in hypothesis-driven studies, underlying the importance of CCASs also in the era of GWAS. In particular, we highlighted the role of rs17849071 and rs17524488, whose association was not significant in previous studies but became statistically significant after increasing the sample size with the present meta-analysis. Overall, our in-depth analysis showed that some a priori hypotheses formulated in previous studies were confirmed and could have realistic bases for shedding some lights in the etiology of DTC. Our GWAS had an adequate statistical power to detect small size effects (>85% of power for SNPs with MAF>0.05, relative risk of 1.4 and type I error α = 0.05), that is reinforced with the data already published through the meta-analysis. However, we cannot exclude that other SNPs could be associated with DTC but failed to replicate in the present study. For example, it is worth mentioning that 13 SNPs (Supplementary Table S3) were found associated with the risk of DTC in the meta-analysis of literature data alone, although all but rs965513 were not significant in the present GWAS. Ideally, all these SNPs should be replicated in a large and independent series of cases and controls to further confirm their involvement in DTC predisposition (24). In conclusion, our findings provide additional evidence that common genetic variants have a role in DTC initiation and/or progression. Further cutting-edge studies, as novel GWASs, next-generation sequencing analysis, fine-mapping or genome-wide interactions studies, are needed to characterize all the predisposing risk factors for DTC.
Disclosure of Potential Conflicts of Interest
No potential conflicts of interest were disclosed.
Authors' Contributions
Conception and design: R. Elisei, C. Romei, K. Hemminki, F. Gemignani, A. Försti, S. Landi
Development of methodology: G. Figlioli, A. Försti
Acquisition of data (provided animals, acquired and managed patients, provided facilities, etc.): G. Figlioli, R. Elisei, F. Bambi, B. Chen, A. Cristaudo, K. Hemminki, A. Försti
Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): G. Figlioli, R. Elisei, O. Melaiu, B. Chen, A. Köhler, K. Hemminki, A. Försti
Writing, review, and/or revision of the manuscript: G. Figlioli, C. Romei, M. Cipollini, A. Cristaudo, F. Gemignani, A. Försti, S. Landi
Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): G. Figlioli, A. Köhler
Study supervision: R. Elisei, K. Hemminki, F. Gemignani, S. Landi
Grant Support
This work was funded by the Istituto Toscano Tumori, grant system 2010.
The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.
Footnotes
Note: Supplementary data for this article are available at Cancer Epidemiology, Biomarkers & Prevention Online (http://cebp.aacrjournals.org/).
- Received June 29, 2015.
- Revision received January 22, 2016.
- Accepted January 23, 2016.
- ©2016 American Association for Cancer Research.
References
Volume 25, Issue 4
