Thyroid function, pernicious anemia and erythropoiesis: a two-sample Mendelian randomization study

Abstract

Autoimmune thyroid disease (AITD) and pernicious anemia (PA) often coexist, but the directionality is unknown. In a two-sample Mendelian randomization (MR) analysis, using summary statistics from large genome-wide association studies (GWASs) in Europeans (N = 49 269–755 406), we examined the genetic associations between thyroid function, PA and markers of erythropoiesis. We performed inverse variance weighted random-effects MR, several sensitivity MR analyses, and bidirectional MR and MR Steiger for directionality. AITD and PA were associated bidirectionally (P ≤ 8 × 10⁻⁶). Neither euthyroid thyroid stimulating hormone (TSH) nor free thyroxine (FT4) were causally associated with PA. One standard deviation (SD) increase in euthyroid FT4 regulated by genetic variants in deiodinases 1 and 2 genes (DIO1/DIO2), corresponding to low-normal free triiodothyronine (FT3) levels, was causally associated with a pernicious/macrocytic anemia pattern, i.e. decreased erythrocyte counts (rank-based inverse normal transformed β = −0,064 [95% confidence interval: −0,085, −0,044], P = 8 × 10⁻¹⁰) and hemoglobin (−0.028 [−0.051, −0.005], P = 0.02) and increased mean corpuscular hemoglobin (0.058 [0.025, 0.091], P = 5 × 10⁻⁴) and mean corpuscular volume levels (0.075 [0.052, 0.098], P = 1 × 10⁻⁸). Meanwhile, subclinical hyperthyroidism mirrored that pattern. AITD was causally associated with increased erythrocyte distribution width (P = 0.007) and decreased reticulocyte counts (P ≤ 0.02), whereas high-normal FT4 regulated by DIO1/DIO2 variants was causally associated with decreased bilirubin (−0.039 (−0.064, −0.013), P = 0.003). In conclusion, the bidirectional association between AITD and PA suggests a shared heritability for these two autoimmune diseases. AITD was causally associated with impaired erythropoiesis and not autoimmune hemolysis. Additionally, in euthyroid individuals, local regulation of thyroid hormones by deiodinases likely plays a role in erythropoiesis.

Introduction

Autoimmune thyroid disease (AITD) is characterized by an autoimmune attack on the thyroid gland with lymphocytic infiltration of the gland and presence of thyroid autoantibodies. AITD often coexists with other autoimmune diseases, like autoimmune gastritis, particularly as part of the autoimmune polyendocrine syndrome type 2 (1,2). In autoimmune gastritis, the gastric epithelial lining is destroyed, leading to vitamin B12 deficiency due to lack of intrinsic factor, and consequently pernicious (macrocytic and megaloblastic) anemia (PA) with decreased erythrocyte counts and hemoglobin levels and increased mean corpuscular volume (MCV) (3–6). In accordance with this, AITD and PA are associated in many clinical and observational studies (7–11). Since AITD and PA are very common, thyroid hormones, hemoglobin and vitamin B12 are among the most frequently ordered blood tests in clinical practice.

In addition to the pathological association between AITD and PA, there is also in vitro and in vivo evidence that thyroid hormones physiologically regulate erythropoiesis (12–14). Bioavailability of thyroid hormones is, among other mechanisms, regulated locally by the deiodinases DIO1 and DIO2, which convert the prohormone thyroxine (T4) into the active hormone triiodothyronine (T3), and DIO3, which converts T4 into inactive reverse T3 (15). Even subclinical and euthyroid thyroid functions are associated with hemoglobin concentrations (16,17). Likewise, treatment of both hypo- and hyper-thyroidism increases hemoglobin levels, suggesting a causal role of thyroid function in anemia (18,19).

It is unknown if a causal directionality between AITD and PA exists, or the association is bidirectional due to an underlying shared heritability independent of thyroid function. The genetic mechanism for the physiological association of thyroid function with erythropoiesis is also unknown. Yet another uncertainty is whether the pathological and the physiological associations of thyroid function with erythropoiesis have a common denominator. An approach to investigate these mechanisms is the Mendelian randomization (MR) design. In a recent MR study, genetically predicted euthyroid TSH and FT4 levels were not associated with anemia defined by decreased hemoglobin levels, except for genetic variants in DIO3OS (DIO3 opposite strand, possibly involved in DIO3 expression levels) gene, suggesting that intracellular regulation of thyroid hormones might play a role in developing anemia (20). However, that study did not investigate AITD, PA or markers of erythropoiesis.

The aim of this study is to investigate: 1) the causal directionality between AITD and PA using bidirectional MR, Steiger filtering and MR Steiger, 2) whether other markers of genetically predicted thyroid function (overt and subclinical hypothyroidism, euthyroid TSH and FT4 levels and subclinical hyperthyroidism) are causally associated with risk of PA, and 3) whether genetically predicted thyroid function is associated with markers of erythropoiesis, i.e. erythrocyte counts and indices, reticulocyte counts and bilirubin. Moreover, we stratify genetic variants associated with FT4 by gene location mapped to DIO1 and DIO2 or DIO3OS.

Results

The study design is illustrated in Figure 1. The main results are summarized in Table 1.

Figure 1

Schematic diagrams illustrating the study design. Mendelian randomization (MR) analyses use genetic variants (Z) as instrumental variables (IVs) to estimate the causal effects of exposures (X) on outcomes (Y). The three core assumptions are shown in (A). The genetic variants used as IVs: 1) must be strongly associated with the exposure, 2) must not be associated with confounders and 3) must not be associated with outcome through any other pathway than the exposure. We used genetic variants associated with thyroid function (exposures) to examine the potential causal effect of thyroid function on pernicious anemia risk and markers of erythropoiesis (outcomes). The summary statistics (β coefficients and standard errors [SEs]) for the IV-exposure (β_ZX) and the IV-outcome (β_ZY) associations were freely available and provided by large GWAS consortia. The overall causal estimates were calculated by performing meta-analyses across the individual causal estimates (Wald ratio: β_XY = β_ZY/β_ZX and their SEs [SE_XY = SE_ZY/β_ZX]). In the bidirectional MR analysis, we used IVs for both autoimmune thyroid disease (AITD) and pernicious anemia in order to evaluate whether AITD causes pernicious anemia, pernicious anemia causes AITD, or whether the associations are bidirectional (B). The dashed boxes represent that AITD and PA are considered as exposures as well as outcomes, depending on the direction of the association examined (arrows).

Table 1

Summary of the main results: causal effects of thyroid function traits on pernicious anemia risk, erythrocyte counts and indices and markers of erythropoiesis

	Thyroid function traits as exposuresa
	Thyroid disease		Euthyroid state
Outcomes	AITDa	Overt hypothyroidisma	High-normal FT4 DIO1/DIO2	Subclinical hyperthyroidism
Pernicious anemia (Fig. 2)	↑↑b	↑b	÷	÷
Erythrocyte counts and indices:	(Supplementary Material, Fig. S1)	(Supplementary Material, Fig. S2)	(Fig. 3)	(Fig. 4)
Erythrocyte count	÷	÷	↓↓	↑
Hemoglobin	÷	÷	↓	÷
MCV	÷	÷	↑↑	↓
Markers of erythropoiesis: (Fig. 5)
Reticulocyte count	↓	↓	÷	÷
Bilirubin	÷	÷	↓↓	÷

↑: increased. ↓: decreased. ÷: no association. Corresponding figures are shown in parentheses.

Single arrow (↓/↑): low power (36–57%) to detect the association. Corresponding P-values for inverse variance weighted random-effects (IVW-RE) Mendelian randomization (MR) varied between 0.003 and 0.03.

Double arrow (↓↓/↑↑): high power (>98%) to detect the association. Corresponding P-values for IVW-RE MR varied between 10⁻⁵ and 10⁻⁹ (except P = 0.01 for FT4 DIO1/DIO2 and bilirubin).

The causal effect sizes are not shown, due to limited clinical use (43).

Power calculations were performed using: https://sb452.shinyapps.io/power/.

^aAITD: autoimmune thyroid disease was mainly overt hypothyroidism, and vice versa.

^bAITD was associated with pernicious anemia (bidirectionally), but not with erythrocyte counts and indices. This suggests that the association between AITD and pernicious anemia is likely due to shared heritability, e.g. autoimmunity, rather than regulation by changes in the hypothalamic-pituitary-thyroid axis per se. AITD was causally associated with impaired erythropoiesis (decreased reticulocyte count) and not autoimmune hemolysis (bilirubin).

High-normal free T4 (FT4, thyroxine) regulated by genes for deiodinases 1 and 2 (DIO1/DIO2), corresponding to low-normal free T3 (FT3, triiodothyronine) was associated with a hematological pattern of pernicious, i.e. macrocytic anemia (↓ erythrocyte count and hemoglobin level, and ↑ MCV). This was mirrored by the opposite hematological pattern in subclinical hyperthyroidism, corresponding to high-normal FT3. Thus, in euthyroid individuals, local regulation of thyroid hormones by deiodinases likely plays a role in erythropoiesis.

Pernicious anemia

Genetic predisposition to AITD was associated with increased risk of PA in the main inverse variance weighted (IVW) analysis (P = 8 × 10⁻⁶, Fig. 2), all sensitivity MR analyses and after Steiger filtering (Supplementary Material, Tables S3 and S4). Between-instrument heterogeneity was high (I² = 58–78%, Supplementary Material, Table S3). There was some evidence of pleiotropy (Egger intercept: P_Egger = 0.003), but low risk of regression dilution bias (I²_GX = 95%, Supplementary Material, Table S3). MR Steiger assessed the direction of the association to be true (Supplementary Material, Tables S3 and S4).

In the bidirectional MR analysis, genetic predisposition to PA was associated with increased risk of AITD in the main IVW (P = 8 × 10⁻⁷) and all sensitivity MR analyses (Supplementary Material, Table S5). Between-instrument heterogeneity was high, but there was no evidence of pleiotropy and low risk of regression dilution bias (I² = 69–97%, P_Egger = 0.98, I²_GX = 97%, Supplementary Material, Table S5). Steiger filtering did not remove any single nucleotide polymorphisms (SNPs), and MR Steiger assessed the direction of the association to be true (Supplementary Material, Table S5).

Genetic predisposition to overt (but not subclinical) hypothyroidism was associated with increased risk of PA in the main IVW (P = 0.03, Fig. 2), Radial MR and MR-PRESSO (MR Pleiotropy RESidual Sum and Outlier) (P < 0.05, Supplementary Material, Table S3), but not in other sensitivity MR analyses nor after Steiger filtering (all P ≥ 0.05, Supplementary Material, Tables S3 and S4). Between-instrument heterogeneity (I²) decreased from 86% to 18% after excluding Radial MR ‘outliers’ (Supplementary Material, Table S3). There was no evidence of pleiotropy and low risk of regression dilution bias (P_Egger = 0.72, I²_GX = 95%, Supplementary Material, Table S3).

Genetically predicted TSH and FT4 (all FT4 SNPs, and FT4 SNPs regulated by DIO1/DIO2 or DIO3OS) and genetic predisposition to subclinical hyperthyroidism were not associated with risk of PA (all P ≥ 0.05, Fig. 2, Supplementary Material, Tables S3 and S4).

Erythrocyte count and indices

Genetic predisposition to AITD was associated with increased erythrocyte (red cell) distribution width (RDW) in the main IVW (P = 0.007) and other sensitivity MR analyses, except MR Egger, but not with erythrocyte counts nor other indices (Supplementary Material, Fig. S1, Tables S3 and S4). Between-instrument heterogeneity was high, but there was no evidence of pleiotropy, and low risk of regression dilution bias (I² = 59–93%, P_Egger = 0.26, I²_GX = 95%, Supplementary Material, Table S3).

Genetic predisposition to hypothyroidism (overt and subclinical) was not associated with erythrocyte counts or indices (Supplementary Material, Figs S2 and S3, Table S3).

One SD increase in genetically predicted TSH was associated with increased mean corpuscular hemoglobin (MCH) (rank-based inverse normal transformed β = 0.026 [95% confidence interval: (0.001, 0.051)], P = 0.04) and MCV (0.026 (0.002, 0.051), P = 0.04) (Supplementary Material, Fig. S4, Table S3). Between-instrument heterogeneity was high, but there was no evidence of pleiotropy or regression dilution bias (I² > 50, P_Egger > 0.73, I²_GX = 98%, Supplementary Material, Table S3).

Genetically predicted FT4 was not associated with erythrocyte counts or indices (Supplementary Material, Fig. S5, Table S3). However, one SD increase in genetically predicted FT4 regulated by DIO1/DIO2 was associated with decreased erythrocyte count (−0,064 [−0,085, −0,044], P = 8 × 10⁻¹⁰) and hemoglobin levels (−0.028 [−0.051, −0.005], P = 0.02) and increased MCH (0.058 (0.025, 0.091), P = 5 × 10⁻⁴) and MCV (0.075 [0.052, 0.098], P = 1 × 10⁻¹⁰) (Fig. 3, Supplementary Material, Table S3). Between-instrument heterogeneity was low-to-moderate (except for MCH). There was no evidence of pleiotropy and low risk of regression dilution bias (I² < 30, P_Egger > 0.08, I²_GX = 96%, Supplementary Material, Table S3). Corresponding causal estimates (β-coefficients) for FT4 regulated by DIO3OS were similar but imprecise (large SDs), and the 95% CIs did not differ from 0 (Supplementary Material, Fig. S6, Table S3).

Genetic predisposition to subclinical hyperthyroidism was associated with increased erythrocyte count (P = 3 × 10⁻³) and decreased MCH (P = 8 × 10⁻⁴) and MCV (P = 1 × 10⁻³) (Fig. 4, Supplementary Material, Table S3). There was no substantial evidence of between-instrument heterogeneity or pleiotropy but a moderate risk of regression dilution bias (I² < 19, P_Egger > 0.36, I²_GX = 59–61%, Supplementary Material, Table S3).

The above results (for erythrocyte counts and indices) were similar in sensitivity MR analyses. Steiger filtering did not remove any SNPs, and the correct causal direction was confirmed by MR Steiger (Supplementary Material, Tables S3 and S4).

Reticulocytes and bilirubin

Genetic predisposition to AITD and overt hypothyroidism was associated with decreased reticulocyte count (P = 0.01 and P = 0.02, respectively, Fig. 5, Supplementary Material, Table S3). There was evidence of moderate-to-high between-instrument heterogeneity as well as pleiotropy, but low risk of regression dilution bias (Supplementary Material, Table S3).

One SD increase in genetically predicted FT4 regulated by DIO1/DIO2 was associated with decreased levels of total and direct bilirubin (−0.039 [−0.064, −0.013], P = 3 × 10⁻³ and − 0.052 [−0.086, −0.017], P = 3 × 10⁻³, respectively, Fig. 5).

The above results (for reticulocytes and bilirubin) were similar in sensitivity MR analyses. Steiger filtering did not remove any SNPs, and the correct causal direction was confirmed by MR Steiger (Supplementary Material, Tables S3 and S4).

No other thyroid function traits were associated with reticulocytes or bilirubin (Supplementary Material, Tables S3 and S4).

Figure 2

Causal effects of thyroid function on risk of pernicious anemia. TSH: thyroid stimulating hormone. FT4: free thyroxine. FT4, DIO1/DIO2: FT4 regulated by genetic variants in deiodinase 1 and 2 genes (DIO1 and DIO2). FT4, DIO3OS: FT4 regulated by genetic variants in deiodinase 3 opposite strand (DIO3OS) gene. AITD: autoimmune thyroid disease. OR: odds ratio. CI: confidence interval.

Figure 3

Causal effects of FT4 regulated by DIO1/DIO2 variants on erythrocyte count and indices. FT4: free thyroxine. DIO1/DIO2: genetic variants in deiodinase 1 and 2 genes. MCHC: mean corpuscular hemoglobin concentration. MCH: mean corpuscular hemoglobin. MCV: mean corpuscular volume. β: beta coefficient. CI: confidence interval. The causal effects per one standard deviation increase in levels of FT4 regulated by DIO1/DIO2 are expressed as rank-based inverse normal transformed units (of erythrocyte count and indices).

Figure 4

Causal effects of subclinical hyperthyroidism on erythrocyte count and indices. MCHC: mean corpuscular hemoglobin concentration. MCH: mean corpuscular hemoglobin. MCV: mean corpuscular volume. β: beta coefficient. CI: confidence interval. The causal effects of subclinical hyperthyroidism, a binary exposure, do not have a clear interpretation, why only the null hypothesis can be tested.

Discussion

This is the first MR study to provide evidence for a bidirectional association between genetic predisposition to AITD and PA, suggesting shared heritability for these two common autoimmune diseases. Genetically high euthyroid TSH, FT4 and FT4 regulated by DIO1/DIO2 or DIO3OS were not associated with PA. However, genetically high euthyroid FT4 regulated by DIO1/DIO2 was directly associated, while genetic predisposition to subclinical hyperthyroidism was inversely associated with a hematological pattern often seen in PA, i.e. decreased erythrocyte count and hemoglobin level, and increased MCH and MCV. Genetic predisposition to AITD was associated with increased RDW and decreased reticulocyte count, whereas genetically high euthyroid FT4, regulated by DIO1/DIO2, was associated with decreased bilirubin levels, suggesting an impaired erythropoiesis.

Conforming to this, an increased RDW (manifested as unequal size of erythrocytes, i.e. anisocytosis) often reflects impaired erythropoiesis and is common in PA and other anemias caused by nutritional deficiencies, i.e. vitamin B12, folate and iron (21).

AITD and PA are both autoimmune diseases (7–10), and in most patients with hypothyroidism or hyperthyroidism, the underlying etiology is attributed to autoimmunity (22,23). The bidirectional association between genetically predicted AITD and PA, as observed in our study, is consistent with an underlying shared heritability between these two diseases.

The fact that neither hypothyroidism (subclinical and overt) nor genetically high euthyroid TSH and FT4 were robustly associated with erythrocyte counts and indices suggests that the association between AITD and PA is not regulated by changes in the hypothalamic-pituitary-thyroid axis per se but rather a shared underlying autoimmunity. In support of this finding, antibodies against gastric parietal cells (24) and intrinsic factor (25), a protein produced by parietal cells and responsible for vitamin B12 absorption, are found in patients with AITD. Thus, an underlying autoimmunity most likely affects both the thyroid (leading to AITD) and autoimmune gastritis, which leads to a cascade of events from lack of intrinsic factor and thereby vitamin B12 deficiency to impaired erythropoiesis and PA. In accordance with this AITD was also associated with increased RDW and decreased reticulocyte counts reflecting impaired erythropoiesis.

We found that genetically high euthyroid FT4, regulated by DIO1/DIO2, was associated with a hematological pattern of pernicious, i.e. macrocytic anemia. This was mirrored by the opposite hematological pattern in subclinical hyperthyroidism. Collectively, these findings suggest that in euthyroid individuals, regulation of erythropoiesis is not controlled centrally (by the hypothalamic–pituitary–thyroid axis) but locally, most likely in the bone marrow, by the balance between FT4 and free triiodothyronine (FT3) regulated by the activity of DIO1 and DIO2 (determined by the genetic variants in the DIO1 and DIO2 genes). FT4 is converted by deiodinases DIO1 and DIO2 to (active) FT3 and converted by DIO3 to (inactive) reverse T3 (15). In erythrocytes, T3 activates the thyroid hormone receptor TRβ, which together with nuclear receptor coactivator NCOA4 regulates the final stages of erythrocyte maturation (12). In our study, the genetic variants in DIO1/DIO2 favored high euthyroid FT4 levels; thus, if the genetic variants were flipped to favor high FT3 levels, the hematological pattern would result in increased erythrocyte count and hemoglobin and decreased MCV, like the pattern we observed for subclinical hyperthyroidism. Thus, the association between FT4 regulated by DIO1/DIO2 variants and erythrocyte counts and indices may just reflect normal up- and down-regulation in the erythropoiesis depending on the body’s needs. We speculate, without any solid evidence, that such a regulation could take place 1) at the entrance to or in the tricarboxylic acid pathway which, as an intermediary step, also provides one of the major substrates for heme synthesis, namely succinyl coenzyme A, and which is also regulated by thyroid hormones (26), and/or 2) by regulating the other major substrate for heme metabolism, namely glycine, which is associated with FT4 (27). However, these speculative mechanisms need to be clarified in further studies.

A severe complication of PA is hemolytic anemia with elevated bilirubin either due to destruction of erythrocytes in the bone marrow or in the circulation (28,29). In our study, genetic predisposition to AITD and overt hypothyroidism were associated with decreased reticulocyte counts, whereas genetically high euthyroid FT4, regulated by DIO1/DIO2, was associated with decreased bilirubin levels, both of which argue against hemolysis. Thus, the most likely explanation for this finding is a decreased synthesis of erythrocytes and hemoglobin, as bilirubin is the terminal product of the heme degradation.

While a recent MR study reported that high euthyroid FT4 regulated by DIO3OS, but not by DIO1/DIO2 variants, was causally associated with increased risk of anemia (defined as decreased hemoglobin levels) (20), we found the opposite, i.e. high euthyroid FT4 regulated by DIO1/DIO2, but not DIO3OS variants, was associated with (a hematological pattern of pernicious) anemia (decreased hemoglobin and erythrocyte counts and increased MCV). This may represent a chance finding in the recent MR study. Alternatively, our study may lack power since the causal estimates for DIO3OS were similar to those for DIO1/DIO2 variants but imprecise as there were only two SNPs mapped to the DIO3OS gene, explaining 0.2% of the variation (r²_exposure) in FT4 levels.

Although genetically high euthyroid FT4 regulated by DIO1/DIO2 and genetic predisposition to subclinical hyperthyroidism were associated with erythrocyte counts and indices, we were unable to examine the association between overt hyperthyroidism and markers of erythropoiesis due to the lack of GWAS variants. Importantly, we would expect the association with erythropoiesis to be in the opposite directions for subclinical and overt forms of hyperthyroidism, since the overt form is observationally associated with anemia (16,18), while we found the subclinical form to be causally associated with increased erythropoiesis. Another potential limitation of our study is that while we decreased the risk of population stratification by restricting the analyses to Europeans only, our findings may not be applicable to other ancestries. Yet another potential limitation that could have biased our estimates is the substantial overlap between the participants in the GWASs included in the two-sample MR, particularly for the AITD and PA, both of which included the UK Biobank (30). However, it is unlikely that this bias would be substantial (30), and the two-sample MR methods can be used even within a single large population sample, such as the UK Biobank (31).

In conclusion, the bidirectional association between AITD and PA suggests a shared heritability for these two common autoimmune diseases. AITD was associated with impaired erythropoiesis but not autoimmune hemolysis. Additionally, in euthyroid individuals, local regulation of thyroid hormones by deiodinases likely plays a causal role in erythropoiesis, i.e. active thyroid hormone FT3 increases erythropoiesis.

Materials and Methods

We used genetic variants associated with thyroid function as instrumental variables (IVs) to examine the potential causal effect of thyroid function on PA risk and markers of erythropoiesis in population-based studies. The conceptual illustration of the MR approach with the three core assumptions is shown in Figure 1A. In the bidirectional MR approach, we used genetic variants associated with both AITD and PA in order to evaluate whether AITD causes PA, PA causes AITD, or whether the associations are bidirectional (Fig. 1B).

Figure 5

Causal effects of thyroid function on reticulocyte counts and total bilirubin levels. FT4: free thyroxine. FT4, DIO1/DIO2: FT4 regulated by genetic variants in deiodinase 1 and 2 genes (DIO1 and DIO2). AITD: autoimmune thyroid disease. β: beta coefficient. CI: confidence interval. The causal effects of AITD and hypothyroidism, i.e. binary exposures, do not have a clear interpretation, why only the null hypothesis can be tested. The causal effects per one standard deviation increase in euthyroid FT4 levels regulated by DIO1/DIO2 are expressed as rank-based inverse normal transformed units (of respective outcomes).

Identification of genetic variants used as instrumental variables

We searched published GWASs for genetic variants associated with thyroid function and PA and included those that: 1) were reported as top independent genetic variants by the published studies (we did not perform additional pruning with respect to linkage disequilibrium, i.e. no additional r² threshold), 2) were autosomal (not on X-chromosome), 3) were SNPs, 4) were associated with thyroid function at GWAS significance level (P ≤ 5 × 10⁻⁸) and 5) had an F-statistic of at least 10, which is commonly considered as sufficient strength threshold for IVs (32). Studies identifying genetic variants for thyroid function are listed in Supplementary Material, Table S1. All ethics approvals were collected by the relevant GWASs this study is based upon.

We identified 89 SNPs associated with AITD (33), 18 SNPs associated with overt (34) and eight with subclinical hypothyroidism (35), 81 SNPs associated with TSH (35,36), 31 SNPs associated with FT4 (four SNPs were mapped to DIO1, three to DIO2 and two to DIO3OS) (35) and seven SNPs associated with subclinical hyperthyroidism (35).

AITD was defined as Graves’ disease (N = 2400), Hashimoto’s thyroiditis (N = 397) and other non- autoimmune hypothyroidism (N = 27 437) in the GWAS of 30 234 cases and 755 172 controls from Iceland and UK Biobank (33).

SNPs for overt hypothyroidism were from the GWAS based on 23andMe data on 17 558 cases with any one of both subclinical and overt hypothyroidism (including thyroid surgery and medication use) and 117 083 controls of European ancestry (34).

Subclinical hypothyroidism (cases = 3440, controls = 49 983) was defined as FT4 within and TSH above the population-based reference range and did not include participants with overt hypothyroidism (35).

TSH and FT4 SNPs were based on a GWAS of individuals with TSH (N = 54 288) and FT4 (N = 49 269) levels within their cohort-specific reference ranges (35). Measures of TSH and FT4 were normalized to within the population-based reference range and scaled to SD units (35).

Subclinical hyperthyroidism (cases = 1840, controls = 49 983) was defined as FT4 within and TSH below the population-based reference range and did not include participants with overt hyperthyroidism (35).

We identified five top independent SNPs associated with PA in 2166 cases and 659 516 controls from population-based biobanks: Estonian Biobank, UK Biobank and FinnGen study (10).

Data access

All data used in this study are already in public domain and accessed as described below.

Outcomes

The outcomes examined were PA, AITD and markers of erythropoiesis, i.e. erythrocyte count and indices (RDW, hematocrit, hemoglobin, mean corpuscular hemoglobin concentration [MCHC], MCH and MCV), reticulocytes (reticulocyte percentage, reticulocyte count, immature reticulocyte fraction and mean reticulocyte volume) and bilirubin (direct and total).

We extracted summary statistics (log-odds and standard errors [SEs]) for the association between thyroid function SNPs and PA from the GWAS of 2166 PA cases and 659 516 controls from population-based biobanks: Estonian Biobank, UK Biobank and FinnGen study (http://www.geenivaramu.ee/tools/pernicious_anemia_Laisketal2021_sumstats.gz) (10). PA was defined as ICD10 code D51.0 in Estonian Biobank and UK Biobank and as D51 (vitamin B12 deficiency anemia) in FinnGen (10).

Summary statistics (odds ratios [OR] and P-values, from which we calculated log-odds and SEs) for the association between the PA SNPs and AITD were from a GWAS of 30 234 AITD cases and 725 172 controls from Iceland and the UK Biobank (33).

We extracted summary statistic (β-coefficients and SEs) for rank-based inverse normalized (irn) erythrocyte count and indices from the Blood-Cell Consortium (BCX) GWAS of up to 563 946 individuals of European ancestry (UK Biobank and several smaller studies) (http://www.mhi-humangenetics.org/en/resources/) (37).

The reticulocyte and bilirubin data sets are a part of the UK Biobank ‘Hematology Data’ release (https://biobank.ndph.ox.ac.uk/showcase/ukb/docs/haematology.pdf) and ‘Biomarker Project’ release (https://biobank.ndph.ox.ac.uk/showcase/showcase/docs/serum_biochemistry.pdf), respectively. We extracted summary statistics for crude and standardized (i.e. irn) outcomes from the rapid UK Biobank GWAS based on the latest release of genetic data, i.e. GWAS round 2 (August 1, 2018) on up to 361 194 individuals of white-British ancestry, provided by the Neale lab (http://www.nealelab.is/uk-biobank/).

Statistical analysis

Analyses were performed in R version 4.0.4 using TwoSampleMR (38,39), RadialMR (40) and MendelianRandomization (41) packages. The R code used in this study is publicly available at: https://github.com/AlisaDK/AITD_PA_erythropoiesis.

For each SNP, the effect allele was defined as the allele associated with increased level of the relevant exposure. For all outcomes, we extracted summary statistics for the relevant exposure SNPs and aligned the effect to the effect allele (data harmonization). This excluded following inconsistent and palindromic SNPs with effect allele frequencies close to 50%: one TSH SNP (rs2979181), three AITD SNPs (rs1257920, rs2271194, rs4293777) and one FT4 SNP (rs11078333). Thus, as IVs, we used 86 SNPs for AITD, 18 SNPs for overt and eight for subclinical hypothyroidism, 80 for TSH, 30 for FT4 (seven for DIO1/DIO2 and two for DIO3OS) and seven SNPs for subclinical hyperthyroidism. Summary statistics (β-coefficients or log-odds and SEs) for all IV-exposure and IV-outcome associations are shown in Supplementary Material, Table S2. However, not all SNPs were present in all GWASs, and since we did not use proxy SNPs, the numbers of IVs across MR analyses vary slightly.

We use the term ‘genetic predisposition’ to AITD, hypo- and hyperthyroidism and PA, as these are binary exposures, and the term ‘genetically predicted’ TSH and FT4, as these are continuous exposures. Importantly, for binary exposures, we could only test the causal null hypothesis (95% confidence intervals and P-values), because the causal estimates do not have a clear interpretation. Therefore, for binary exposures, only the P-value but not the causal estimates are shown in the Results section (42). Furthermore, even for the continuous exposures, the size of the causal estimates is of limited clinical use (43). The main analysis was the IVW multiplicative random-effects meta-analysis applied across the individual instrumental estimates and their SEs (44). The heterogeneity of the individual causal estimates was assessed by the Cochran’s Q and corresponding I² index (range: 0–100%, high I² means high heterogeneity) (45). We used the random-effects model, since it provides a more conservative estimate (same point estimate, but increasing confidence intervals with increasing heterogeneity) than the fixed-effects model of the meta-analysis. However, both IVW models assume that all SNPs are valid IVs. Furthermore, since the IVW estimates are calculated using the ‘first-order’ IVW, the IVW methods make the assumption that the IV-exposure association is measured without error. A strong violation of the NO Measurement Error (NOME) assumption may induce the regression dilution bias in the IVW estimate towards the null hypothesis (46). We therefore performed sensitivity analyses with different assumptions regarding IV validity and pleiotropy: Radial MR, MR-PRESSO, weighted median (WM) and MR-Egger regression analyses.

In order to minimize horizontal pleiotropy, we excluded outlier SNPs with greatest heterogeneity contributions, by employing the Radial MR (40,46) as well as the comparable simulation-based approach in MR-PRESSO (47). We re-calculated the IVW (random-effects) estimates after removing the ‘outliers’ identified by Radial MR using modified second-order weights at a Bonferroni corrected α = 0.05/N_SNPs (see ‘IVW, RadialMR outliers’ in the Results) (40,46).

The WM estimate is defined as the median, i.e. 50th percentile, of the IVW empirical distribution function of the Wald ratio estimates (48). Thus, WM relaxes the first IV assumption, as it assumes that at least 50% of the weight contributed by genetic variants comes from valid IVs (48).

While all the causal estimates can be obtained as slopes of the linear regression of the IV-outcome on the IV-exposure association estimates, only the MR-Egger regression line is not constrained to zero. Since the MR-Egger intercept (with corresponding P-value for intercept = 0) can be interpreted as the average pleiotropic effect across all IVs, we used it to test for directional pleiotropy (49). Furthermore, we quantified the NOME violation by the I²_GX (range: 0–100%, low I²_GX means high risk of regression dilution bias, i.e. a shift of causal estimates towards the null hypothesis) (49).

The directionality was assessed by the bidirectional MR, as well as by Steiger filtering and MR Steiger. Steiger filtering is based on the assumption that each valid IV should explain more variance in the exposure than the outcome (r²_exposure > r²_outcome) and therefore removes all genetic variants in conflict with this assumption (38,39). Meanwhile, MR Steiger uses all genetic variants and assesses the causal direction to be true if all the genetic variants combined explain more variance in the exposure than the outcome, and false if opposite is the case (38).

We considered P < 0.05 as the significance level threshold. We specifically did not use the Bonferroni correction, as it assumes independence within the exposures, within the outcomes and within the analyses. This is overly conservative due to the high degree of correlation within the thyroid traits (including AITD), within the markers of erythropoiesis (including PA), as well as within the MR analyses performed (all based on meta-analysis). Furthermore, a P-value above a chosen threshold is not a proof of lack of an association (50). We therefore performed an overall evaluation of each result with a P < 0.05 individually by considering the effect size (and SE), biological plausibility and consistency across the examined exposures, outcomes and MR analyses.

Acknowledgements

The authors thank the GWAS consortia: ThyroidOmics, 23andMe, UK Biobank (including Neale lab), Icelandic Biobank, Estonian Biobank, FinnGen and BCX for providing the summary statistics data used in this study.

Conflict of Interest statement. The authors declare no known competing financial interests or personal relationships that could have appeared to influence this study.

Funding

A.D.K. is funded by an unrestricted grant from Novo Nordisk Foundation.