Autoimmune thyroid diseases (AITD), mainly comprising autoimmune thyroiditis/Hashimoto’s disease (HT) and Graves’ disease (GD), are complex polygenic disorders with similarities and differences, primarily programmed by genetic factors(Reference Antonelli, Ferrari and Ragusa1,Reference Zhang, Liu and Li2) . AITD is one of the most common pathologies in endocrine diseases and is the most common cause of hypothyroidism and hyperthyroidism. AITD affects 2–5 % of the population, with a majority being women(Reference McLeod and Cooper3). AITD patients are often accompanied by other autoimmune disorders, one of which is celiac disease (CeD). CeD is an immune-mediated disease caused by ingestion of gluten, affecting 1 % of the total population(Reference Ludvigsson, Bai and Biagi4,Reference Conrad, Misra and Verbakel5) , yet many patients with CeD go undiagnosed(Reference Lebwohl, Sanders and Green6–Reference Lohi, Mustalahti and Kaukinen9). Its pathogenesis involves the immune response to gluten protein, the activation of immune cells and the damage to the intestinal mucosa. A lifelong strict gluten-free diet (GFD) is the key to treatment(Reference Lebwohl, Sanders and Green6).
In recent years, AITD co-morbid with CeD has attracted attention. A systematic review and meta-analysis of twenty-seven studies calculated the median prevalence of CeD in AITD to be 1·4 %, with a higher prevalence of CeD in hyperthyroidism/Graves’ patients (2·6 %)(Reference Roy, Laszkowska and Sundström10,Reference Ferrari, Fallahi and Ruffilli11) ; the incidence of AITD was also found to be increased in CeD patients, especially HT(Reference Dore, Fanciulli and Rouatbi12,Reference Bibbò, Pes and Usai-Satta13) . This co-morbid phenomenon opens up new possibilities for the care of AITD: Is it possible to migrate the GFD treatment to patients with AITD? It seems feasible, as it has been hypothesised that gluten may contribute to the deterioration of immune-mediated disorders(Reference Lerner, Shoenfeld and Matthias14), and GFD treatment exhibits a general anti-inflammatory effect in extraintestinal autoimmune inflammatory diseases(Reference Bruzzese, Scolieri and Pepe15,Reference Passali, Josefsen and Frederiksen16) . Therefore, GFD is very popular in AITD, especially HT, despite the lack of favourable evidence to suggest its benefit in AITD(Reference Pobłocki, Pańka and Szczuko17,Reference Malandrini, Trimboli and Guzzaloni18) .
Mendelian randomisation (MR) analysis serves as a viable substitute for randomised controlled trials in determining causality. This methodology, grounded in the principles of Gregor Mendel’s laws of segregation and independent assortment, is employed to emulate a trial that is randomised on a genetic level. Researchers utilise genetic variants, which have a strong association with the exposure under investigation, as instrumental variables (IV). This enables the estimation of the causal effect of a particular exposure or intervention on a specific outcome while reducing the impact of confounding variables. Moreover, the application of MR allows to effectively bypass the problem of reverse causation bias, considering that genetic variations are determined.
In this study, we aim to explore the causal relationship between CeD and AITD, particularly with HT, using the method of MR. Additionally, we attempt to provide genetic evidence that GFD is beneficial for AITD.
Methods
Data source
Genetic instruments were retrieved from the MRC IEU OpenGWAS database (https://gwas.mrcieu.ac.uk) and the NHGRI-EBI GWAS catalogue (https://www.ebi.ac.uk/gwas/), which are the two largest databases of genome-wide association study (GWAS) analysis summary statistics, as well as from previously published GWAS studies. GWAS we used were corrected for age, age2, sex, age*sex, age2 * sex and the first twenty principal components of genotype data. All the data we used were selected from publicly available databases and do not require ethical approval. GWAS data sources for each exposure are provided in Table 1.
Overview of the data sources
Autoimmune thyroid diseases
Summary GWAS data for HT and GD were sourced from a study by Sakaue et al. (Reference Sakaue, Kanai and Tanigawa19), which investigated 220 human phenotypes, including AITD. This study integrated UK Biobank (UKB) and FinnGen data and included 395 640 Europeans. The UKB project recruited approximately 500 000 people aged between 40 and 69 years across the UK between 2006 and 2010. This study retrieved clinical information and then defined HT and GD using phecodes. The FINNGEN project had the participation of six regional and three country-wide Finnish biobanks to track the health status of participants from birth to death information by linking to the national health registries national health registry (1969–2016); the study used pooled data from FinnGen version 3 data and defined diseases with ICD-10 standards. We retained the GWAS data from the European population.
Celiac disease and gluten-free diet
CeD and GFD data are from the study by Jiang, L. et al.(Reference Jiang, Zheng and Fang20), where the original data were derived from the UKB project. For CeD, the researcher used UKB’s unsorted text data (UKB datafield 41 202). UKB defines CeD using the ICD-10 code K90·0, and the researcher then grouped the ICD-10 code into Phecode based on the ICD-10 to Phecode v.1·2 map20. The final CeD data presented meets the PheCode 557·1 standards definitions. For GFD, the UKB project collected data from 58 913 volunteers from 2009 to 2012 who responded to a 24-h recall diet questionnaire, ‘Types of Special Diets Followed,’ via online follow-up (More specific information is on https://biobank.ctsu.ox.ac.uk/crystal/field.cgi?id=20086 ).
Assumptions
This study adheres to the three assumptions of MR research: (1) Genetic variations are closely associated with the risk factor (relevance). (2) Genetic variations are independent of confounding factors (independence). (3) Genetic variations influence the outcome only through the risk factor (exclusivity)(Reference Emdin, Khera and Kathiresan21). A schematic diagram is shown in Fig. 1.
A schematic diagram of the three main assumptions in this Mendelian randomisation analysis.
Method for selecting instrumental variables
First, meaningful SNP were selected from the GWAS summary data used as exposure (selection criteria: for AITD, P < 5E-08, linkage disequilibrium r 2 = 0·001, and the width of the linkage disequilibrium region is 10 000 kb, to ensure that each SNP is independent and to exclude the impact of gene pleiotropy on the results; for CeD and GFD, to obtain sufficient data for analysis, we expanded the selection criteria: P < 5E-06, linkage disequilibrium r 2 = 0·01, and the width of the linkage disequilibrium region is 5000 kb). We calculated the F-statistics to quantify the strength of genetic variation, with F < 10 being a weak IV that needs to be excluded. The formula for calculating the F-statistics is F = (N-K-1)/K × R2/(1-R2), where N is the sample size of the exposure database, K is the number of SNP and R 2 is the proportion of variation explained by the SNP in the exposure database. The formula for calculating R 2 is R 2 = 2 × MAF × (1 – MAF) × β 2, where MAF is the minor allele frequency and β is the allele effect value. Then, the exposure-related SNP screened out from the GWAS summary data of the study outcome were extracted; missing SNP were directly excluded without using alternative methods. SNP directly related to the study outcome were excluded (P < 5E-08). Before the MR analysis, we performed an MR-PRESSO analysis to exclude any outliers with potential pleiotropy to ensure the reliability of our MR estimates(Reference Verbanck, Chen and Neale22). According to the above steps, the remaining SNP were finally used as genetic instruments.
Methods of Mendelian randomisation analysis
The primary analysis method is the inverse variance weighted (IVW) method(Reference Burgess, Bowden and Fall23), which combines Wald ratio estimates to obtain a consistent estimate of the causal impact of exposure on the outcome, characterised by not considering the existence of an intercept term during regression. We also used other MR analysis methods, including MR Egger regression(Reference Slob, Groenen and Thurik24), weighted median estimator(Reference Bowden, Davey Smith and Haycock25), simple mode(Reference Hartwig, Davey Smith and Bowden26) and weighted mode(Reference Hartwig, Davey Smith and Bowden26) to verify the causal relationship between exposure and outcome. The main difference between the MR Egger regression method and the IVW analysis method is the addition of an intercept term, mainly used to judge whether there is horizontal pleiotropy. The weighted median estimator assumes that at least half of the IV in the analysis are valid and then obtains a consistent estimate of the causal relationship. Simple mode classifies SNP according to causal effects and divides similar values into a cluster, using the cluster with the most SNP to evaluate the estimated causal effect. Weighted mode weights the causal effect value of SNP in each cluster with the number of SNP in the module, and the returned result is a temporary estimate with the maximum SNP number weight. The use of a consistent causal effect estimate of the weighted mode requires the ‘ZEMPA assumption’ (zero modal pleiotropy assumption)(Reference Hartwig, Davey Smith and Bowden26).
Sensitivity analysis
In this study, Cochran’s Q test was used to judge the heterogeneity of the IVW model; when P < 0·05, heterogeneity is considered to exist, and at this time, the random effects model of IVW is used for causal inference. MR Egger regression analysis was used to detect whether genetic pleiotropy exists; when the intercept term differs significantly from 0, it indicates the presence of horizontal pleiotropy. Leave-one-out sensitivity test was performed to check deeper.
Software and package
Statistical analysis was conducted using the two-sample MR package (version 0.5.7) in R software (version 4.3.0).
Result
Instrumental variables
The IV were selected according to the method previously described. Details of IV for exposures are shown in online Supplementary Table S1.
Mendelian randomisation with gluten-free diet as exposure
We first explored the GFD as exposure in relation to CeD and AITD. Complete primary and sensitivity analyses of MR are presented in Fig. 2 and Table 2.
Mendelian randomisation for estimating the association between gluten-free diet, celiac disease and autoimmune thyroid disease using inverse variance weighted method, weighted median estimator, MR Egger regression, simple mode and weighted mode.
Sensitivity analyses of Mendelian randomisation for each outcome in relation to gluten-free diet
There was no evidence to support a genetic liability related to a GFD with CeD (ORIVW = 1·178 (95 % CI 0·983, 1·412), P = 0·077), HT (ORIVW = 0·996 (95 % CI 0·936, 1·059), P = 0·895) and GD (ORIVW = 1·007 (95 % CI 0·976, 1·039), P = 0·648) (Fig. 2). There was no indications of heterogeneity in the MR Egger and IVW regression tests, with no evidence of any horizontal pleiotropy in the MR Egger intercept (Table 2).
Mendelian randomisation with celiac disease as exposure
Complete primary and sensitivity analyses of MR for each exposure in relation to AITD are presented in Fig. 3 and Table 3.
Mendelian randomisation for estimating the association between celiac disease and autoimmune thyroid disease using inverse variance weighted method, weighted median estimator, MR Egger regression, simple mode and weighted mode.
Sensitivity analyses of Mendelian randomisation for celiac disease as exposure
In the aforementioned MR analysis, there was no statistically significant effect of CeD (ORIVW = 1·007 (95 % CI 0·976, 1·039), P = 0·648) and genetic susceptibility associated with HT (Fig. 3). Under the MR Egger intercept, the results were also robust, but they did not meet the criteria of the heterogeneity test (Table 3). Given that we used random effects IVW as the primary MR analysis, heterogeneity in this context is not unacceptable(Reference Burgess, Davey Smith and Davies27).
Our results support a genetic liability related to CeD (ORIVW = 1·060 (95 % CI 1·013, 1·109), P = 0·012) with GD (Fig. 3), and this result was robust under the MR Egger intercept, but they did not meet the criteria of the MR Egger and IVW heterogeneity test. Therefore, we repeated the MR analysis with a more stringent selection of CeD-related IV (selection criteria: P < 5E-06, LD clumping window = 10000 kb, LD clumping r 2 = 0·001, effective IV are shown in the online Supplementary Table S1) and obtained negative results. These results passed the MR Egger, IVW and MR Egger intercept tests (Table 4).
Mendelian randomisation (MR) and sensitivity analyses of MR for celiac disease in relation to Graves’ disease with a more stringent selection
Mendelian randomisation with autoimmune thyroid disorder as exposures
Complete primary and sensitivity analyses of MR for reverse causation between CeD and AITD are presented in Fig. 4 and Table 5.
Mendelian randomisation for estimating the association of reverse causation between celiac disease and autoimmune thyroid disease using inverse variance weighted method, weighted median estimator, MR Egger regression, simple mode and weighted mode.
Sensitivity analyses of Mendelian randomisation for reverse causation between celiac disease and autoimmune thyroid disease
In the reverse MR, strong evidence suggests that HT plays a pathogenic role in increasing the risk of CeD (ORIVW = 1·544 (95 % CI 1·153, 2·068), P = 0·004) (Figure 5 shows the scatter plot). There is no evidence of heterogeneity in our MR Egger (P = 0·064) and IVW regression tests (P = 0·097), nor is there any evidence of horizontal pleiotropy in the MR Egger intercept (P = 0·729) (Table 5). The results remain robust under further leave-one-out analysis (Fig. 6).
Scatter plot for Hashimoto’s disease on celiac disease.
Leave-one-out for Hashimoto’s disease on celiac disease.
However, there is no evidence to suggest that GD has a pathogenic or protective effect (ORIVW = 1·189 (95 % CI 0·988, 1·428), P = 0·066) on CeD (Fig. 4), and under the MR Egger intercept, this result is robust (Table 5).
Discussion
Our study provides the following evidence: (1) An increased risk of HT is associated with an increased risk of CeD, and this causal relationship is unidirectional, that is, suffering from HT increases the risk of CeD, whereas there is no adequate evidence suggesting that developing CeD increases the risk of HT. (2) An increased risk of GD is not associated with an increased risk of CeD. (3) There is no evidence to support that a GFD is beneficial in reducing the risk of AITD.
Both CeD and AITD have dysregulated cellular immunity, and the leaky gut state triggered by celiac disease permits the passage of non-self epitopes and antigens, which may trigger autoimmunity(Reference Lerner, Jeremias and Matthias28). It is in light of this significant correlation that GFD to improve AITD has become increasingly popular, but whether or not patients benefit is controversial. Previous studies have found potential beneficial effects of GFD in AITD: In female HT patients with IgA-tTG seropositivity but asymptomatic CeD, GFD could decrease the gradient of thyroid peroxidase antibody and thyroglobulin antibody(Reference Krysiak, Szkróbka and Okopień29); thyroid function abnormalities were reversed in some newly diagnosed CeD patients after 1 year of GFD(Reference Valentino, Savastano and Tommaselli30). The results of a meta-analysis also support the beneficial effects of a GFD(Reference Piticchio, Frasca and Malandrino31). However, some studies have provided contrary evidence, with no effect on the levels of TPO antibodies in ten cases (37 %) of newly diagnosed CeD patients after 1 year of GFD treatment, and even with the establishment of GFD in patients without CeD, the thyroid volume still significantly decreased(Reference Metso, Hyytiä-Ilmonen and Kaukinen32). These conflicting opinions drew our attention. In clinical practice, the strict implementation of GFD can be brutal for clinical patients, especially Asian patients, and the lack of evidence can again shake clinical patients’ resolve to implement GFD; therefore, obtaining firm evidence is critical. Considering conducting clinical trials of GFD is difficult, especially with low patient compliance and results confounded by multiple factors, we used an MR method, ‘nature’s randomised trial’, trying to answer the question: Is there any relationship between GFD and CeD and AITD? Our hypothesis was that ‘GFD improves AITD by improving CeD’, but surprisingly, our results did not support this idea, while supporting the conclusions of the clinical study by Metso S et al.(Reference Metso, Hyytiä-Ilmonen and Kaukinen32).
However, the negative result is still instructive: although CeD and AITD share cellular immune dysregulation, they may not interact, but rather, more likely, there is a more sophisticated network of regulatory mechanisms, that is, pathogenic and protective roles in parallel, resulting in a failure to show a clear pathogenic or protective role at the disease phenotype level. This suggests that further studies at the cellular as well as molecular level are necessary and will be the next step in our in-depth research, and the positive finding of an increased risk of CeD due to HT will be an important entry point. Besides, HT and GD, as two subtypes of AITD, are potentially coexisting or transforming into each other. In the clinic, it is able to detect positivity of relevant antibodies and changes in different antibody titres before and after transformation, indicating that the immune mechanisms that play a major role in the disease network also change during the transformation process. Although both belong to AITD, our results only clarify the possible pathogenic role of HT in CeD; this result suggests that the abnormal cellular immune state in HT is more relevant to CeD. There is a greater possibility to uncover the core pathways and targets from this positive result.
Our study possesses multiple strengths. Given that lifestyle elements and socio-economic status are susceptible to real-world bias, it becomes challenging to accurately determine their actual impact through observational studies. The utilisation of MR analysis minimises the potential influence of residual confounding factors. In addition, considering the heterogeneity of HT and GD, our MR analysis independently analyses these two AITD subtypes, making it one of the few known MR analyses that independently analyse HT. The data we used are based on the currently available largest GWAS summary statistics for HT. To enhance the robustness of the results, thorough sensitivity analyses were performed using multiple MR methods, and no violations of any assumptions were found. Practically, our study results help provide more dietary strategy information for AITD patients.
This investigation also presents certain limitations. Initially, we confined our research participants to individuals of European ancestry in an effort to mitigate population stratification bias, and thus, the findings are not yet generalisable to other ethnic populations. We attempted to extrapolate the study results to the Asian population but failed due to the inability to find GWAS data related to CeD and GFD. We look forward to more GWAS studies on the Asian population. Second, in the choice of dietary exposure, we also tried other GWAS data, but it was difficult to continue the analysis due to the lack of sufficient valid IV. Third, similar to all MR studies, we are unable to tackle the issue of unobserved horizontal pleiotropy. Fourth, MR analyses, the ‘natural RCT’, use genetic variants strongly associated with exposure factors as IV to assess the causal relationship between exposure factors and outcomes, thus avoiding interference from common confounding factors such as acquired environmental and social factors and excluding reverse causal effects. Therefore, the order of diagnosis does not influence our results, nor does age or gender, but we cannot exclude that age and gender play some mediating role. Fifth, the original data for CeD were obtained from the UKB project, where the majority of cases were identified by serology and endoscopy, both of which are widely recognised as reliable diagnostic criteria for CeD, although a proportion of cases were based on symptomatic presentation alone, which may lack the specificity and accuracy of serological or histological confirmation (More specific information is on https://biobank.ndph.ox.ac.uk/showcase/field.cgi?id=21069 ). Differences in diagnostic methods, such as relying on symptoms alone v. strict serologic or histologic confirmation, can lead to heterogeneity in case definition, which creates limitations when interpreting and comparing results. Sixth, the inclusion of GFD was largely based on a single self-reported piece of data, which may be subject to recall bias, social desirability bias or individual differences in understanding. And the questionnaire only addresses the past 24 h of diet, making it difficult to assess whether individuals have adhered to a GFD over time. It is difficult to carry out a strict assessment on a large scale, which is precisely why we chose MR analyses to assess the causal relationship between exposure factors and outcomes, by using genetic variants closely related to the exposure factors as IV, thus avoiding the interference of common confounding factors. However, we also expect more objective and standardised methods, such as biomarker-based assessments or validated dietary questionnaires, to provide a more rigorous raw database for analysis.
The findings of this study support that an increased risk of HT is associated with an increased risk of CeD but do not support a protective effect of a GFD on HT. This suggests that the ‘effect’ of GFD in HT may be due to the improvement of CeD-related symptoms in patients with HT. However, some studies have shown that adhering to a GFD may positively impact the absorption of selenium and vitamin D elements necessary for thyroid function and health(Reference Stazi and Trinti33–Reference Verma, Lata and Khanna35). We cannot rule out that GFD may be beneficial to patients with HT through these effective intermediary factors but is masked by other non-beneficial factors. In the future, continuous large-scale GWAS summary statistics of dietary pattern genetic variations are needed to more accurately assess the effective dietary strategies for AITD.
Supplementary material
For supplementary material/s referred to in this article, please visit https://doi.org/10.1017/S0007114525103917.
Acknowledgements
We thank the associate editor and the reviewers for their useful feedback that improved this paper.
This study was supported by the National Natural Science Foundation of China (82274486).
P-Y. H.: Data curation, writing – original draft. L-C. Y.: Data curation, writing – review and editing. J-Y. Y.: Data curation, writing – review and editing. L-Z. H.: Data curation, writing – review and editing. Z-B. X.: Writing – review and editing, supervision. Y-R. S.: Supervision, writing – review and editing.
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
