Unique genetic and risk-factor profiles in clusters of major depressive disorder-related multimorbidity trajectories

Abstract

The heterogeneity and complexity of symptom presentation, comorbidities and genetic factors pose challenges to the identification of biological mechanisms underlying complex diseases. Current approaches used to identify biological subtypes of major depressive disorder (MDD) mainly focus on clinical characteristics that cannot be linked to specific biological models. Here, we examined multimorbidities to identify MDD subtypes with distinct genetic and non-genetic factors. We leveraged dynamic Bayesian network approaches to determine a minimal set of multimorbidities relevant to MDD and identified seven clusters of disease-burden trajectories throughout the lifespan among 1.2 million participants from cohorts in the UK, Finland, and Spain. The clusters had clear protective- and risk-factor profiles as well as age-specific clinical courses mainly driven by inflammatory processes, and a comprehensive map of heritability and genetic correlations among these clusters was revealed. Our results can guide the development of personalized treatments for MDD based on the unique genetic, clinical and non-genetic risk-factor profiles of patients.

Fig. 1. Overview of the main methods.

A Rationale and hypothesis of the study. Accumulating evidence suggested that MDD is frequently comorbid not only with other psychiatric disorders but also with several somatic diseases contributing to worse health-related outcomes and decreasing quality of life^–. Thanks to network medicine and system biology approaches, it has been demonstrated that comorbid conditions partially represent common biological mechanisms^–. Furthermore, directly related comorbidities of depression, where the relationships are not mediated by other disorders, represent stronger molecular-level relationships and are time-dependent (i.e., vary with onset age). Finally, a recent comorbidity mapping study of asthma supported that comorbidities are indeed suitable to delineate distinct subgroups of complex multifactorial disorders. B The cohort-specific datasets contain the onset ages of diseases in three-character ICD-10 categories. Data were collected from the participants over various periods, depicted by the length of the grey lines, with disease onsets marked by an ‘x’. Participant trajectories were discretized into cumulative time intervals, as shown at the bottom of the figure. C The structure of the inhomogeneous dynamic Bayesian network used. The boxes correspond to intervals, the nodes in the boxes correspond to diseases, and the solid and dashed edges indicate direct relations between the diseases. This method determined the strongly relevant MDD-related multimorbidities; these nodes are in the Markov boundary of the target variable, indicated by the grey-shaded region and a thick black node border. Genetic and other non-genetic variables also influenced the onset of the diseases (dotted edges). One aim of the study was to identify pleiotropic genetic variants (edges with α) that influence the onset of MDD and its related multimorbidities. These variants confound the direct relationship (edge β) between MDD and its strongly relevant comorbid conditions. D Overview of the study pipeline. We determined MDD-related cross-cohort clusters of all participants in the UKB, CHSS, and THL cohorts by utilizing the temporal trajectories of the participants’ MDD-related multimorbidity burden. The seven identified clusters were then characterized based on disease and non-genetic risk-factor profiles and genetic contributions, and the findings were validated in the two independent cohorts (the FinnGen and SHIP cohorts).

Fig. 2. Temporal disease patterns in the clusters according to UKB (N = 502,504).

A The average onset ages of the most prevalent (>5%) cross-cohort diseases (per line) according to the seven clusters in the UKB cohort. The node colour indicates the clusters, and the node size is proportional to the observed prevalence of the disease (in %) in the cluster. MDD is highlighted with red colour. B The onset distribution of all cross-cohort diseases in the clusters of the UKB cohort is shown as violin plots, where the colour indicates the cluster. Each of the seven violin plots represents one cluster, computed using data from the entire cohort (N = 502,504), with individuals weighted by their posterior probability of belonging to the respective cluster. The onset distribution and the mean age of MDD onset are indicated with inline red violin plots and red vertical bars, respectively.

Fig. 3. Disease risk in the clusters of the discovery and validation cohorts.

A Coefficient of cluster membership (hazard ratio [HR] in the weighted Cox proportional hazards regression model) with respect to the onset of each cross-cohort disease in the UKB cohort (N = 502,504). The top five diseases with the strongest increase/decrease in risk in each cluster are indicated in the plot and listed on the right. The colour of the markers corresponds to the main ICD-10 disease category. D: Diseases of the blood and blood-forming organs, E: Endocrine, nutritional and metabolic diseases, F: Mental, behavioural and neurodevelopmental disorders, G: Diseases of the nervous system, H: Diseases of the eye, ear and mastoid process, I: Diseases of the circulatory system, J: Diseases of the respiratory system, K: Diseases of the digestive system, L: Diseases of the skin and subcutaneous tissue, M: Diseases of the musculoskeletal system and connective tissue, N: Diseases of the genitourinary system. B Values and 95% confidence intervals of cluster membership coefficients (hazard ratios) from the weighted Cox proportional hazards regression models for the onset of MDD across various cohorts. Points indicate coefficient values and error bars represent the 95% confidence intervals. Colours represent the different cohorts. C Weighted Kaplan‒Meier estimates of MDD-free survival in the various cohorts throughout participants’ lifespans. Survival curves are labelled by cluster numbers, and the colours of the curves indicate the distinct clusters. The dotted grey curves indicate the mean MDD-free survival in the whole cohort, regardless of cluster membership. (In A and C: UKB N = 502,504; THL N = 41,092; CHSS N = 645,913; FinnGen N = 385,640; SHIP N = 1449).

Fig. 4. Results from GWAS analyses in the UKB (N = 249,167) and FinnGen cohorts (N = 277,252).

A Gene-based genome-wide Manhattan plots for the seven clusters in the UKB cohort. Association analyses were first performed for each cluster using linear regression to test the association between each SNP and the posterior log odds of cluster membership, controlling for age, sex, the first ten genetic principal components, and the genotyping array. Next, MAGMA gene-level analysis was performed to identify putative significant genes using a SNPwise-multi model, defining the SNP set of each gene with a ± 10 kb window. In the plot, nominal p-values are displayed. The genome-wide significant genes are indicated with red dots, and the significance threshold (2.7 × 10^-6) is depicted with a dashed dark red line. B Genetic correlation (r_g) plot from GWAS summary statistics on the posterior log-odds of cluster membership among Clusters 1–7 in the UKB cohort. The colour of the dots indicates the value of the genetic correlation. C Genetic correlation (r_g) plot from GWAS summary statistics on the posterior log-odds of cluster membership among Clusters 1–7 in the FinnGen cohort. The area and the colour of the circles represent the magnitude and direction (blue = positive, red = negative) of the genetic correlation between two clusters. D Overlap between genome-wide significant genes from MAGMA analyses of Clusters 1–7 in the UKB cohort. Black dots indicate clusters within the comparison. The intersection size corresponds to the number of genes uniquely shared by these clusters.

Fig. 5. Non-genetic risk-factor profiles for each cluster in the UKB (N = 249,167) and THL cohorts (N = 23,786).

Simple linear regression models, including one factor at a time with age and sex as covariates, were calculated for each cluster in the A UKB cohort and B THL cohorts. The posterior log-odds of being in a given cluster were the dependent variable. The direction of the triangles reflects the sign of the coefficient (upwards = positive; downwards = negative), and the colour reflects the magnitude. Statistical analyses were performed using two-sided t-tests to assess the significance of each factor’s effect on cluster membership. Adjustments for multiple comparisons were made using the Bonferroni correction. The size of the triangles is proportional to the −log₁₀ p-value, and only significant values are shown (−log₁₀(p) > 4). Sex was coded as follows: 1 = male, 2 = female. The risk factors of stress and neuroticism score were not available in the THL cohorts. *Alcohol intake and depression score were not available in the FinHealth17 and Finrisk cohorts, respectively, but were available in the other THL cohort.