Bidirectional Mendelian Randomization and Multiphenotype GWAS Show Causality and Shared Pathophysiology Between Depression and Type 2 Diabetes

Type 2 diabetes is a disease characterized by chronic hyperglycemia, while depression is its frequent comorbidity, potentially because of shared risk factors (1). It has been shown that depression, even at subclinical levels, increases the risk of incident type 2 diabetes by 25–60% (2), whereas others have shown that type 2 diabetes increases the risk of depression by 40–60% (3).

The causality of the associations from observational studies remains unclear because of unmeasured confounding and potential reverse causation. However, this could be circumvented in part through Mendelian randomization (MR), an approach that assesses potential causality between phenotypes using genetic variants as instruments, since genes are allocated randomly at birth and are free from confounding (4). To date, only one MR study from China reported a possible causal link between type 2 diabetes and depression (5). Other studies have demonstrated a causal effect of depression on type 2 diabetes and the components of the metabolic syndrome (6,7). However, the reported causal estimates in one of these studies may be biased because of the use of nonoverlapping samples (8) in its MR studies. Additionally, although Tang et al. (6) aimed to address the role of BMI in the causal relationship between depression and type 2 diabetes by using BMI-adjusted summary statistics, they did not formally estimate the indirect effects mediated via BMI.

Previous large-scale genome-wide association studies (GWAS) for type 2 diabetes and depression reported 403 and 102 associated genomic loci for these diseases, respectively (9,10). However, standard GWAS investigate each disease independently, without considering the genetic correlation between related phenotypes and their heritabilities.

Here, we addressed the causal relationship between depression and type 2 diabetes by conducting an MR study using summary statistics from nonoverlapping samples of previous GWAS of depression (10) and type 2 diabetes (9). Additionally, we used the UK Biobank (UKBB) to perform a multiphenotype GWAS (MP-GWAS) for type 2 diabetes and depression to identify shared genetic loci. For depression, we compared two assessment approaches, clinically diagnosed major depressive disorder (MDD) and depressive symptoms based on self-report.

To assess causality between type 2 diabetes and depression, we performed a two-sample bidirectional MR, first using depression as a risk factor and type 2 diabetes as an outcome, then the reverse, from two nonoverlapping data sets. We then used a two-step MR approach to assess the potential mediation via BMI in the depression and type 2 diabetes relationship (11) (Supplementary Fig. 1). The single nucleotide polymorphisms (SNPs) used as genetic instruments for type 2 diabetes, depression, and BMI were from previous large-scale European GWAS meta-analyses (9,10,12). To limit potential bias introduced by sample overlap, we used the type 2 diabetes summary statistics without the UKBB data set since UKBB was part of the depression meta-analyses.

All MR analyses were conducted using the R software package TwoSampleMR version 0.5.4 (13). In the type 2 diabetes GWAS summary statistics used, 96 SNPs (inclusive of 6 proxy variants with a minimum r² ≥0.8) out of the 102 independent (r² <0.01) depression SNPs were available. We excluded seven palindromic SNPs (A/T or C/G) with intermediate allele frequencies (minor allele frequency >45%) to ensure that the effects of the SNPs for the two phenotypes were aligned to the same forward strand allele. The genetic instruments for type 2 diabetes included 403 genetic SNPs associated from a previous GWAS (9). In the depression summary statistics, 343 (inclusive of 4 proxy variants with a minimum r² ≥0.8) out of 403 independent (r² <0.01) type 2 diabetes SNPs were available for the analysis. To obtain the total causal effect estimate, we applied the inverse variance weighted (IVW) method (4). We performed sensitivity MR analysis using weighted median (WM) and MR-Egger regression methods to evaluate the potential violations of the MR assumptions (14) (Supplementary Methods). F-statistic was used to evaluate the instrument strength, where F >10 indicates the presence of a strong instrument. The F-statistic indicated a good instrument strength for both type 2 diabetes (F = 32.87) and depression (F = 42.99). We assessed heterogeneity between the causal estimates from each SNP using Cochran Q test. The sensitivity of causal inference to any individual genetic variant was tested by leave-one-out analysis. We used the Strengthening the Reporting of Observational Studies in Epidemiology Using MR (STROBE-MR) reporting guideline to facilitate readers’ evaluation of our results (15) (Supplementary Table 1).

To assess the potential indirect effect via BMI in the depression-to-type 2 diabetes total effect, we additionally performed a two-step MR (11) using 97 BMI SNPs (12). In the first step, the causal effect of depression (exposure) on BMI (mediator) was determined. In the second step, we estimated the causal effect of BMI on type 2 diabetes (outcome). The indirect effect via BMI as a mediator was then obtained by taking the product of the effects from steps 1 and 2 of the two-step MR (Supplementary Methods). The causal estimates for both steps were obtained using the IVW method, while the WM was used as a sensitivity analysis.

We used data from the UKBB (https://www.ukbiobank.ac.uk), which consists of 500,000 individuals aged between 40 and 69 years at recruitment. Genetic data were available for 488,377 individuals in the UKBB genotyped using the UKBB BiLEVE array (n = 49,979) and the UKBB Axiom Array (n = 438,398) (16). The genetic data were imputed using the Haplotype Reference Consortium, the UK10K, and phase 3 of the 1000 Genomes Project (16), resulting in ∼90 million variants available for association testing.

The FinnGen data set summary statistics were used for replication of our GWAS and MP-GWAS results (https://www.finngen.fi/fi). The June 2022 data freeze used in our analysis comprised 309,154 individuals. Summary statistics for 3,095 phenotypes are publicly available for analysis. FinnGen study participants were genotyped using the Illumina and Affymetrix chip arrays (Illumina Inc. [San Diego, CA] and Thermo Fisher Scientific [Santa Clara, CA, https://www.thermofisher.com). The data were then imputed using the Sequencing Initiative Suomi version 3 imputation panel (https://sisuproject.fi), resulting in 16,962,023 variants available for association analysis.

In the UKBB, individuals who self-reported to have a diabetes diagnosis by a physician, did not commence insulin medication within 1 year after diagnosis, and were at least 40 years old at the time of diagnosis were classified as type 2 diabetes case subjects. Individuals not meeting these criteria were classified as control subjects. For both case and control subjects, we excluded individuals with gestational diabetes mellitus (field 4041, code 1), those <40 years old at the time of diagnosis (field 2976) and those on insulin medication within the first year of diagnosis (field 2986). Sex-discordant individuals (genotype vs. reported sex) were also excluded from the analyses. We included individuals of all ancestries in the entire UKBB to maximize on sample size. In total, 19,344 case subjects and 463,641 control subjects were available (Supplementary Table 2). The GWAS summary statistics in the FinnGen study were based on the ICD-10 code E11 for type 2 diabetes in 304,769 individuals (49,303 case subjects, 255,466 control subjects).

We defined depression in two ways: ICD-10–coded MDD based on linked data from hospital records, and self-reported depressive symptoms using the Patient Health Questionnaire 9 (PHQ-9) (17).

Individuals with a primary diagnosis of a depressive episode (ICD-10 code F32) and recurrent depression (ICD-10 code F33) were defined as MDD case subjects (hereafter referred to as MDD). Individuals who answered no to the questions, “Have you ever seen a general practitioner for nerves, anxiety, tension or depression?” (field 2090) and “Have you ever seen a psychiatrist for nerves, anxiety, tension or depression?” and no to either “depressed/down for a whole week” (field 4598) or “ever unenthusiastic/disinterested for a whole week” were defined as control subjects. Participants were excluded from the study if they had a diagnosis of bipolar disorder (ICD-10 codes F30 and F31), mixed and other personality disorder (ICD-10 code F61), or schizophrenia (ICD-10 code F20). Participants taking any of 58 antipsychotic medications (field 20003) were also excluded. In total, 5,262 case subjects and 86,275 control subjects were used for the MDD phenotype (Supplementary Table 1). The proportion of individuals with MDD who had a type 2 diabetes diagnosis are shown in Supplementary Table 3.

FinnGen study data had ICD-10 codes F32 and F33 available for MDD. The GWAS summary statistics on MDD were based on a total of 305,192 individuals (33,812 case subjects, 271,380 control subjects).

In UKBB, self-reported depressive symptoms over the previous 2 weeks (from the time of study enrollment) were assessed using the PHQ-9 (17) (Supplementary Table 4). It has been shown that the PHQ-9 is invariant between people with and without diabetes (18), suggesting that its interpretation is similar for both diabetes case subjects and control subjects. Individuals missing responses for more than three PHQ-9 items were excluded from the analysis. Missing PHQ-9 responses for the remaining individuals were imputed using imputeSCOPA software (https://github.com/ImperialStatGen/imputeSCOPA), which implements a random forest approach to impute missing items. The variables sex, age, education qualification, BMI, Townsend Deprivation Index (an area-based measure of deprivation), genotyping array, and eight principal components were included in the imputation model to improve the predictive accuracy of the imputation (Supplementary Table 5). The sum of all nine PHQ-9 items after imputation for all individuals was used for quantitative association analysis. PHQ-9 data were available for 153,079 individuals (Supplementary Table 1). The proportion of individuals with PHQ-9 data and a type 2 diabetes diagnosis are shown in Supplementary Tables 3 and 6. Symptom-based depression phenotypes were unavailable in the FinnGen replication data set.

We performed separate GWAS for type 2 diabetes, MDD, and PHQ-9 in the UKBB data with BOLT-LMM using a linear mixed model (19). We adjusted for age, sex, array, BMI, and the first eight principal components. We analyzed common variants (minor allele frequency >5%), with imputation scores >0.4, Hardy-Weinberg equilibrium P > 1 × 10⁻⁶ and per-SNP variant missingness <0.015. All analyses were performed on human genome build 37. The statistical threshold for genome-wide significant SNPs (signals) used was P < 5 × 10⁻⁸.

We used multitrait analysis of GWAS (MTAG) (20), which implements a generalized IVW- meta-analysis, to increase the power for locus identification, improve SNP effect size estimates for type 2 diabetes and depressive phenotypes, and identify potential multiphenotype genetic variant effects. We used summary statistics of the individual GWAS as inputs of MTAG. In UKBB, two MTAG models were tested, one with type 2 diabetes and MDD and the other with type 2 diabetes and total PHQ-9 scores, to assess the consequence of using two different depression definition criteria (disease diagnosis vs. disease symptoms) in an MP-GWAS approach. To assess the robustness of the phenotype-specific MP-GWAS results, MTAG computes a maximum false discovery rate (maxFDR) statistic, a theoretical upper bound limit on the FDR for a GWAS (20). A maxFDR <5% indicates robust results.

To explore and identify target genes of the identified loci, we used several expression quantitative trait locus (eQTL) databases. We extracted eQTL data from the Genotype-Tissue Expression (GTEx) Portal (https://gtexportal.org) for SNPs identified in our MP-GWAS, focusing on type 2 diabetes and depression relevant tissues (i.e., brain, muscle, liver). In addition, as GTEx does not include data from pancreatic islets, we used recent eQTL data from the TIGER Data Portal (https://tiger.bsc.es) obtained from >500 brain-dead organ donor islets (21). We extracted data for the seven shared SNPs identified in our MP-GWAS from both eQTL studies using a nominal significance threshold of P < 0.05.

Furthermore, we also used GTEx version 7 transcriptome data (22) from European individuals to identify eQTLs using our MP-GWAS summary statistics. We focused on relevant tissues in 1) type 2 diabetes, namely liver, whole pancreas, muscle, subcutaneous adipose, adrenal gland, and whole blood, and 2) depression, including putamen basal ganglia, hippocampus, substantia nigra, frontal cortex, amygdala, and anterior cingulate cortex. For each tissue, the predicted expression levels were then correlated with type 2 diabetes and PHQ-9 MTAG summary statistics. P values were corrected for multiple testing using Bonferroni correction based on the number of genes tested per tissue (Supplementary Table 7). Genes where <80% of the SNPs used in the model were found in the GWAS summary statistics were excluded because of low reliability of association results. This analysis focused on type 2 diabetes and PHQ-9 phenotypes only, as the MDD phenotype was underpowered in both GWAS and MP-GWAS.

Our MR analysis revealed that depression was causally and positively associated with type 2 diabetes using the IVW method, with an odds ratio (OR) of 1.26 (95% CI 1.11–1.44, P = 5.46 × 10⁻⁴). This result was consistent with the WM sensitivity analyses, which showed an OR of 1.18 (95% CI 1.02–1.38, P = 0.0304) (Fig. 1, Supplementary Fig. 2A, and Supplementary Table 8). The MR-Egger test showed no evidence of directional pleiotropy (P = 0.15), further confirming the validity of the results. The Cochran Q statistic for heterogeneity was significant for the IVW method (Q = 158.18, P = 1.02 × 10⁻⁷). However, the leave-one-out analysis showed no extreme outliers, suggesting that the observed association was not changed after removing any single variant (Supplementary Fig. 3). We found no evidence of causality in the reverse direction between type 2 diabetes and depression in the primary or sensitivity analysis (IVW: OR 0.99 [95% CI 0.99–1.01], P = 0.92) (Supplementary Fig. 2B and Supplementary Table 9).

We observed evidence of an indirect effect via BMI (OR 1.09 [95% CI 1.000083–1.19], P = 0.049) using the IVW method, suggesting potential mediation of the depression to type 2 diabetes causality via BMI. The indirect effect was consistent using the WM method (OR 1.09 [95% CI 1.03–1.14], P = 0.0021). The estimated proportion of the mediated effect via BMI was 0.37 (95% CI 0.12–0.58), indicating that any potential mediation via BMI would account for 37% of the total effect. Mediation via BMI was not explored from type 2 diabetes to depression, as we observed no significant total causal effect in this direction.

To identify whether type 2 diabetes and depression have a shared genetic etiology, we first performed a GWAS for both phenotypes separately using the UKBB. For type 2 diabetes, we identified 92 independent SNPs at 84 loci, of which 82 were replicated in the FinnGen data set with nominal significance (P < 0.05) and consistent in direction of effect. For depression, we found three independent SNPs and loci for PHQ-9 and no SNPs associated with the binary depression MDD trait in the UKBB. None of the GWAS SNPs identified in type 2 diabetes were shared with depression (Supplementary Fig. 4 and Supplementary Tables 10–12).

To improve the power to identify shared genetic variants between the two phenotypes, we performed the largest MP-GWAS to date of type 2 diabetes and depressive phenotypes in the UKBB using the respective GWAS summary statistics. The MP-GWAS model with type 2 diabetes and MDD did not identify any significant associations for the binary definition of depression (MDD). For type 2 diabetes, we identified 71 independent signals at 66 loci (Supplementary Figs. 5 and 6 and Supplementary Table 10). The maxFDR for the type 2 diabetes and MDD MP-GWAS results was 1.5% and 7.7%, respectively (Supplementary Table 13), indicating that the MDD results were likely inflated by the higher powered type 2 diabetes GWAS. In FinnGen, the maxFDRs were 11.5% and 25.5%, respectively, indicating highly inflated results for both traits.

In contrast, in the MP-GWAS model with type 2 diabetes and PHQ-9, we identified eight independent SNPs for PHQ-9 (compared with only three SNPs detected in GWAS) (Fig. 2 and Supplementary Table 14), suggesting greater power for SNP discovery. Only the CACNA2D2 locus was reported in both GWAS and MP-GWAS analyses for PHQ-9, although with different, yet highly correlated lead SNPs (CACNA2D2, chromosome 2, SNP_GWAS rs35335661, SNP_MP-GWAS rs1467916, r² = 0.99). For type 2 diabetes, MP-GWAS identified 53 SNPs at 50 loci (compared with 92 identified in the single-phenotype GWAS), of which 24 loci were previously reported (10). We replicated 43 type 2 diabetes signals in the FinnGen cohort (Supplementary Table 14).

In total, we found seven SNPs shared between type 2 diabetes and PHQ-9: rs7903146 (TCF7L2), rs7766070 (CDKAL1), rs1359790 (SPRY2), rs16860235 (IGF2BP2), rs76895963 (CCND2-AS1), rs2972144 (IRS1), and rs10811662 (CDKN2B-AS1) (Table 1 and Fig. 2). The maxFDR for type 2 diabetes and PHQ-9 after MP-GWAS were 0.98% and 1.8%, respectively (Supplementary Table 13), indicating that the results are robust and not influenced by the sample sizes of either GWAS.

To explore whether the seven SNPs shared between type 2 diabetes and self-reported depression had a downstream functional impact, we first extracted eQTL data of the seven identified SNPs using the 1) GTEx Portal in relevant tissues, including muscle, liver, and brain, and 2) the TIGER Portal for pancreatic islets. We found that of the seven shared SNPs between type 2 diabetes and depression, six were associated with the expression of nearby genes in relevant tissues (Table 2). This includes the rs2972144-G risk allele associated with a decreased expression of IRS1 and RP11-395N3.2 in visceral and subcutaneous adipose tissue and RP11-395N3.1 in subcutaneous adipose tissue. Additionally, the rs76895963-T risk allele was associated with decreased expression of CCND2 in pancreatic islets, brain cerebellum, skeletal muscle, and subcutaneous adipose tissue and decreased expression of its antisense CCND2-AS1 in pancreatic islets, brain cerebellum, basal ganglia, and cortex in addition to CCND2-AS2 in the cerebellum. For pancreatic islets, we found an increased expression of NDUFA9 (Table 2).

In addition, using our MP-GWAS summary statistics, we tested in GTEx whether type 2 diabetes and PHQ-9 were associated with shared gene expression changes in several target tissues, except for pancreatic islets as GTEx does not include these data. This analysis identified additional target genes associated with both type 2 diabetes and PHQ-9 in eight tissues, consistent in direction of effect for both phenotypes: subcutaneous adipose (IRS1, NCR3LG1, RP11-395N3.2), amygdala (HSPA1B), frontal cortex (CDKAL1), hypothalamus (EIF2S2P3), skeletal muscle (HLA-DRA, RP11-370C19.2), substantia nigra (BETL1), and whole blood (EIF2S2P3, HLA-DRB1) (Supplementary Fig. 7 and Supplementary Tables 7 and 15–17). Altogether, our data demonstrate a functional impact of our loci in several related tissues.

We performed a large, comprehensive study to investigate the relationship between type 2 diabetes and depression and found evidence for a causal positive association from depression to type 2 diabetes, with evidence of mediation by BMI. We also performed an MP-GWAS of the two diseases, highlighting seven shared loci that target nearby genes in several target tissues.

Our MR analyses provide evidence for causality on the previously reported epidemiological associations from depression to type 2 diabetes (2) but not in the reverse direction. These findings are consistent with the pathophysiology for these diseases, whereby depression starts in adolescence or early adulthood (23), whereas type 2 diabetes usually develops later (24), and poor health habits among individuals with depression, including smoking, physical inactivity, and increased caloric intake (and associated overweight), that are known to facilitate the development of type 2 diabetes (2,25). Indeed, our results indicate mediation by BMI in the causal relationship between depression and type 2 diabetes, further supporting this hypothesis. Furthermore, antidepressants frequently induce weight gain, leading to type 2 diabetes (25), and the systemic inflammation associated with increased stress hormone levels, such as cortisol in the context of depression, also favors insulin resistance (26).

While epidemiological studies also show increased risk of depression in people with type 2 diabetes (3), we did not find evidence for causality in this relationship. Our evaluation of instrument strength for the MR analysis suggests that the lack of this causal association is not simply an issue of statistical power. A likely hypothesis is that the epidemiologically observed association is confounded via other factors, which are not easily assessed in epidemiological studies, such as psychosocial factors related to the management of treatment-related pain, i.e., from routine self-injections with insulin, in this middle-age chronic disease. Diabetes distress, the psychological burden caused by dealing with having diabetes and having to care for it, has been linked to depression (27). However, a large international longitudinal study of 14 countries showed that only depressive symptoms, rather than MDD, were predicted by diabetes distress 1 year after diagnosis (28).

Our MP-GWAS revealed seven shared SNPs between type 2 diabetes and self-reported depressive symptoms consistent in their directions of effect. These SNPs mapped to loci implicated in insulin secretion pathways, consistent with previous studies (29,30). We did not find any shared SNPs associated with the binary MDD clinical diagnosis and type 2 diabetes, which we speculate could be due to PHQ-9 being a continuous definition of depression and likely improves power compared with the binary MDD definition. Previous studies have also shown that depression defined based on self-reported symptoms (minimal phenotyping) rather than strict diagnostic criteria enables greater power for locus discovery in GWAS (31).

Our eQTL analysis provided clues about the underlying mechanisms linking the two diseases. For instance, we found that the rs76895963-T risk allele was associated with the decreased expression of CCND2 in the brain, pancreatic islets, and insulin target tissues (adipose and skeletal muscle), suggesting a pleiotropic effect of this locus. CCND2 encodes cyclin D2, which is involved in pancreatic β-cell proliferation and insulin secretion (32). Mouse knockout models of CCND2 have no brain neurogenesis and showed mild depression-like symptoms (33). Therefore, CCND2 could have a multisystem and pleiotropic effect that potentially mediates the relationship between type 2 diabetes and depression.

In addition, we found that the rs2972144-G risk allele was associated with decreased expression of IRS1, which encodes insulin receptor substrate 1, in adipose tissue. IRS1 is a key signaling molecule associated with insulin resistance (34), a condition linked to both type 2 diabetes and depression (35).

Lower expression of CDKAL1 in the frontal cortex was associated with both type 2 diabetes and depressive symptoms. Variation in CDKAL1 has been associated with depressive phenotypes (29). Although we did not find an eQTL for this locus in pancreatic islets, variation at CDKAL1 has been implicated in type 2 diabetes and shown to reduce insulin secretion (36).

We also found two target genes implicating the HLA region (HLA-DRA and HLA-DRB1) in blood and skeletal muscle. The role of depression treatments that target the immune system is very active (37). Type 2 diabetes is also known to have an impact on the immune system, with high blood glucose levels causing an inflammatory response (38).

On the basis of our results, we draw some clinical significance. Antidepressant treatment offered to people with depression at risk for type 2 diabetes should favor those that provide better glycemic control, such as selective serotonin reuptake inhibitors (39). Additionally, people with depression should be encouraged, as part of routine clinical care, to engage in positive lifestyle habits, such as increased physical activity, adequate sleep, and a proper dietary regimen.

This study exhibits several strengths. The two-sample bidirectional MR incorporated nonoverlapping data sets to limit bias and was performed using the MP-GWAS for the two diseases. Using two-step MR, we assessed the role of BMI as a potential mediator in the depression-to-type 2 diabetes causal effect, an approach that has not been used in previous studies. However, there are some limitations to be considered. Our analyses involving MDD in MR, GWAS, and MP-GWAS require validation in larger MDD data sets as they become publicly available. Our study used Hospital Episode Statistics data available in the UKBB to define the diagnosis of MDD (ICD-10 codes F32 and F33). Although not as accurate as gold-standard diagnostic interviews, Hospital Episode Statistics data can, with reasonable certainty, ascertain a diagnosis of unipolar depression (40). Combining the diagnostic categories of F32 (major depressive episode) and F33 (recurrent major depression) may have further increased the heterogeneity of depression in our cohort. To effectively probe the credibility of the MTAG MP-GWAS results, replication using another MP-GWAS method is needed. Furthermore, we could only replicate the type 2 diabetes-MDD analyses because of unavailability of publicly available depressive symptoms GWAS in the FinnGen cohort.

In conclusion, our results highlight the importance of an efficient prevention of type 2 diabetes at the onset of depressive symptoms. Maintaining a healthy weight in the context of its effect on depression and type 2 diabetes comorbidity is also highly recommended.

Acknowledgments. Summary statistics for type 2 diabetes were made available by the Diabetes Genetics Replication and Meta-analysis (DIAGRAM) Consortium. Summary statistics of depression were made available by the Psychiatric Genomics Consortium. The eQTL data analyzed in this study were obtained through the TIGER Data Portal by the T2D Systems Consortium and the GTEx portal. The authors thank all the investigators for sharing the data used in this study. This research was conducted using the UKBB resource under application no. 35327. The authors thank the participants and staff of the UKBB for their valuable contributions, as well as the participants and investigators of FinnGen study. The graphical abstract was designed using images from https://www.flaticon.com.

Funding. This research was funded in part by the Diabetes UK (British Diabetic Association no. 20/0006307), European Union’s Horizon 2020 research and innovation program (H2020 Science with and for Society) (LONGITOOLS, H2020-SC1-2019-874739), Agence Nationale de la Recherche (PreciDIAB, ANR-18-IBHU-0001), and European Union through the Fonds européen de développement regional (FEDER), Conseil Régional des Hauts-de-France (Hauts-de-France Regional Council), and Métropole Européenne de Lille (MEL, European Metropolis of Lille).

Duality of Interest. No potential conflicts of interest relevant to this article were reported.

Author Contributions. J.G.M. wrote the first draft of the manuscript. J.G.M., Z.B., A.N., A.B., P.F., A.K., I.Pr., and M.K. revised the manuscript. J.G.M., Z.B., I.Pu., A.U., M.B., A.K., and M.K. conducted the analyses, provided analytical support, and generated/collected data. Z.B., A.N., P.F., A.K., I.Pr., and M.K. designed the study. All authors helped to interpret the data and read the final version of the manuscript and approved its submission. M.K. is the guarantor of this work and, as such, had full access to all the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.