The co-occurrence of multiple chronic diseases, known as multimorbidity [1], affects 1 in 3 adults. Its prevalence rises with age, affecting 60% of individuals aged between 65 and 74 years and escalating to 80% among those aged ≥85 years [2]. Due to its association with poor prognosis, functional impairment, and reduced quality of life, multimorbidity is considered a global health care challenge [3,4] tied to complex clinical situations, leading to increased encounters with health care professionals, hospitalizations, and pharmacological prescriptions, resulting in a substantial rise in health care costs [5]. The emergence of multimorbidity is not arbitrary and frequently aligns with shared risk factors and underlying pathophysiological mechanisms [6-8] that result from complex interactions between genetic and environmental factors throughout the life span [9]. Perceiving diseases not in isolation but as integral components of a more extensive, interconnected system within the human body has led to the emergence of network medicine [10,11]. Network medicine analyzes disease co-occurrence patterns, aiming to understand the complex connections between diseases to uncover biomarkers, therapeutic targets, and potential interventions [12,13]. The studies investigating the temporal patterns of disease concurrence, or disease trajectories [14,15], rely on a pragmatic approach to this concept to yield a better understanding of the time-dependent relationships among diseases and establish a promising landscape to identify disease-disease causal relationships.

According to this paradigm, a disease-centered approach might lead to suboptimal treatment of patients with multiple chronic conditions, triggering the need to implement new tools to enhance the effectiveness of health services [16]. In this regard, multimorbidity-adjusted health risk assessment (HRA) tools [17-21], such as the morbidity groupers, are crucial for assessing the comprehensive health needs of patients with multimorbidity [22]. HRA uses algorithms and patient data to categorize individuals by risk, aiding health care professionals in customizing interventions, optimizing resource allocation, and enhancing patient outcomes through preventive care. HRA tools facilitate efficient case-finding and screening processes [23]. Case finding targets the individuals at high risk, which is crucial for specialized health care programs, whereas patient screening detects latent illnesses early, enabling cost-effective interventions to prevent disease progression and reduce health care demands.

However, despite their widespread use, prevailing population-based HRA tools such as the Adjusted Clinical Groups [24], Clinical Risk Groups [25], or Adjusted Morbidity Groups (AMG) [4,21] still do not incorporate information on disease trajectories in their calculations. The AMG system is currently used in Catalonia (Spain; 7 million inhabitants) for health policy and clinical purposes. Adding disease-disease association information to the AMG (or other morbidity groupers) may open new avenues for implementing epidemiological impact analyses concerning clusters of interrelated diseases. In addition, it may facilitate the construction of risk groups that accurately represent probabilities of developing new associated disease conditions [26] susceptible to early prevention.

While acknowledging current limitations, this study sought to explore the feasibility of incorporating procedures relevant to the study of disease trajectories [14,15] and novel techniques for analyzing dependency relationships between concomitant diseases [27,28] to improve the capabilities of the current morbidity groupers. This approach might better adjust the estimations of the burden of morbidity to clusters of diseases and improve the ability to anticipate the progression of multimorbidity.

We used major depressive disorder (MDD; F32-F33 in the International Classification of Diseases, 10th Revision, Clinical Modification [ICD-10-CM] [29]) as a use case due to its clinical relevance in multimorbidity management. However, this study pursued to showcase a methodology applicable beyond MDD, allowing for the assessment of the impact of multimorbidity across different clusters of diseases.

This paper describes the process of development and assessment of the Multimorbidity-Adjusted Disability Score (MADS) through an observational retrospective multicentric cohort study, showcasing a pioneering approach that integrates advanced techniques for analyzing disease associations, insights from the analysis of disease trajectories, and a comprehensive scoring method aimed at evaluating the disease burden. The MADS was designed to stratify patients with different health needs according to (1) the disease burden caused by MDD and its comorbidities on individuals and health systems and (2) the risk of morbidity progression and the onset of MDD comorbid conditions.

On the basis of the temporal disease maps among MDD and highly prevalent disease conditions [30] generated using Bayesian direct multimorbidity maps (BDMMs) [27,28], a promising method for filtering indirect disease associations in the context of the European Research Area on Personalized Medicine project “Temporal disease map based stratification of depression-related multimorbidities: towards quantitative investigations of patient trajectories and predictions of multi-target drug candidates” (TRAJECTOME) [31], we combined the probabilities of relevance (PRs) among MDD and its comorbid conditions with the disability weights (DWs) [32], documented in the 2019 revision of the Global Burden of Disease (GBD) study, to compute the MADS. We used the MADS to generate a risk pyramid and stratify the study population into 5 risk groups using different percentiles of MADS distribution. Finally, we analyzed the correspondence between the MADS risk groups and health outcomes through a cross-sectional analysis of mortality and use of health care resources and a longitudinal analysis of disease prevalence and incidence of new disease onsets. The clinical relevance of the identified risk groups was assessed through a multicentric assessment of the findings. To this end, MADS performance was analyzed using data from 3 independent European cohorts from the United Kingdom, Finland, and Spain including >1 million individuals.

The development and evaluation of the MADS involved the following steps (Figure 1).

Step 1 involved computing age-dependent disease-disease PRs using the BDMM method in 4 age intervals (0-20 years, 0-40 years, 0-60 years, and 0-70 years). This analysis resulted in an inhomogeneous dynamic Bayesian network that determined the PR for MDD against the most prevalent co-occurring diseases in the 3 European cohorts considered in TRAJECTOME, namely, the Catalan Health Surveillance System (CHSS) [33], the UK Biobank (UKB) [34], and the Finnish Institute for Health and Welfare (THL) [35] cohorts. The THL cohort amalgamates information from the FINRISK [36] 1992, 1997, 2002, 2007, and 2012; the FinHealth [37] 2017; and the Health [38] 2000 and 2011 studies.

In step 2, combining the PR of every disease condition assessed in the study with its corresponding DW extracted from the GBD 2019 study, we estimated the morbidity burden caused by MDD and its comorbid conditions. The MADS was computed following a multiplicative combination of the PR and DW of all the disease conditions present in an individual.

Step 3 involved using the MADS to stratify patients into different risk levels corresponding to different percentiles of the population-based risk pyramid of each patient cohort: (1) very low risk (percentile ≤50), (2) low risk (percentile 50 to percentile 80), (3) moderate risk (percentile 80 to percentile 90), (4) high risk (percentile 90 to percentile 95), and (5) very high risk (percentile >99).

Finally, in step 4, the correspondence between the MADS risk strata and health outcomes was analyzed through a cross-sectional analysis of use of health care resources, mortality, pharmacological burden, and health care expenditure and a longitudinal analysis of disease prevalence and incidence of new disease onsets. The results were validated through a multicentric replication of the findings in the 3 study cohorts including 1,041,014 individuals.

BDMMs were used to assess direct and indirect associations between MDD and 86 potential comorbid conditions. The set of 86 disease conditions considered in the study had a prevalence of >1% in all the study cohorts. The list of diseases and their associated ICD-10-CM [29] codes are shown in Multimedia Appendix 1.

This step considered information on disease diagnosis (disease conditions were cataloged using the first 3 characters of the ICD-10-CM codes), age at disease onset (the age at disease onset corresponds to the first diagnosis in a lifetime for each ICD-10-CM code), sex, and socioeconomic status (annual average total household income [before tax with copayment exemption] as a categorical variable with 3 categories: <€18,000 [US $19,565.30], €18,000-100,000 [US $19,565.30-108,696], and >€100,000 [US $108,696]).

BDMM analysis resulted in an inhomogeneous dynamic Bayesian network, which was used to compute temporal PR, ranging from 0 (no association) to 1 (strong association), for MDD in conjunction with sex, socioeconomic status, and the set of 86 predetermined diseases [30]. To construct the trajectories, the PR was calculated in 4 different age ranges: 0 to 20 years, 0 to 40 years, 0 to 60 years, and 0 to 70 years. The PRs calculated and used for MADS computation are reported in Multimedia Appendix 1. Further details regarding the core analysis conducted in TRAJECTOME can be found in the study by Juhasz et al [30].

The MADS was developed by weighting the DWs of single diseases according to their estimated PRs against MDD. DWs indicate the degree of health loss based on several health outcomes and are used to indicate the total disability caused by a certain health condition or disease. Often, the DWs present specific disability scores tailored to the severity of the disease. The disease categories, severity distribution, and associated DWs used in this study were extracted from the GBD 2019 study and are reported in Multimedia Appendix 1.

DWs were extracted and aggregated as follows. First, we considered only the DWs of MDD and the set of 86 disease codes. Second, we considered the DWs of all the chronic conditions diagnosed in patients’ lifetime, whereas as the disability caused by acute illnesses is transitory, the DWs for the acute diseases diagnosed >12 months before the MADS assessment were arbitrarily set to 0 (no disability). Third, due to the unavailability of information on the severity of diagnoses, we determined the DWs of each disease condition by calculating the weighted mean of the DWs associated with the disease severity categories and their prevalence. In instances in which the severity distribution was not available, we computed the arithmetic mean of the DWs of each severity category. Fourth, we weighted the DWs according to the PR of each disease condition with respect to MDD. The PRs were adjusted according to the age of disease onset, discretized in the following intervals: 0 to 20 years, 20 to 40 years, 40 to 60 years, and >60 years.

As the DWs do not account for multimorbidity in their estimates, the use of DW independently can cause inaccuracies in the burden of disease estimations, particularly in aging populations that include large proportions of persons with ≥2 disabling disease conditions [39]. Consequently, we combined the DW and the PR for all the disease conditions present in 1 individual following a multiplicative approach (equation 1) [40], aggregating several DWs in a single score that accounts for the overall disability caused by numerous concurrent chronic conditions in which every comorbid disease increases the utility loss of a patient, although it is less than the sum of the utility loss of both diseases independently.

Once calculated, the MADS was used to stratify patients in different levels of risk according to the percentiles of its distribution in the source population, producing the following risk pyramid: (1) very low risk (percentile ≤50), (2) low risk (percentile 50-percentile 80), (3) moderate risk (percentile 80-percentile 90), (4) high risk (percentile 90-percentile 95), and (5) very high risk (percentile >99).

The clinical relevance of the risk strata was assessed through two interconnected analyses: (1) a cross-sectional analysis of health outcomes and (2) a longitudinal analysis of disease prevalence and incidence of new onsets.

The study was conducted using data from 3 public health cohorts.

As a multicentric study, TRAJECTOME accessed data from multiple cohorts, all subject to the legal regulations of their respective regions of origin, and obtained the necessary approvals from the corresponding ethics committees.

For the CHSS cohort, the Ethical Committee for Human Research at Hospital Clínic de Barcelona approved the core study of TRAJECTOME on March 24, 2021 (HCB/2020/1051), and subsequently approved the analysis for the generation and the assessment of the MADS on July 25, 2022 (HCB/2022/0720).

The UKB cohort received ethics approval from the National Research Ethics Service Committee North West–Haydock (reference 11/NW/0382).

The THL cohort integrates information from the FINRISK databases (1997 [ethical committee of the National Public Health Institute; statement 38/96; October 30, 1996], 2002 [Helsinki University Hospital, ethical committee of epidemiology and public health; statement 87/2001; reference 558/E3/2001; December 19, 2001], 2007 [Helsinki University Hospital, coordinating ethics committee; Dnro HUS 229/EO/2006; June 20, 2006], and 2012 [Helsinki University Hospital, coordinating ethics committee; Dnro HUS 162/13/03/11; December 1, 2011]), the FinHealth 2017 (Helsinki University Hospital, coordinating ethics committee; 37/13/03/00/2016; March 22, 2016), and the Health 2000 to 2011 databases (ethical committee of the National Public Health Institute, 8/99/12; Helsinki University Hospital, ethical committee of epidemiology and public health, 407/E3/2000; May 31, 2000, and June 17, 2011).

The ethics committees exempted the requirement to obtain informed consent for the analysis and publication of retrospectively acquired and fully anonymized data in the context of this noninterventional study.

All the data were handled in compliance with the General Data Protection Regulation 2016/679, which safeguards data protection and privacy for all individuals in the European Union (EU). The study was conducted in conformity with the Helsinki Declaration (Stronghold Version, Brazil, October 2013) and in accordance with the protocol and the relevant legal requirements (Law 14/2007 on Biomedical Research of July 3).

The results of the cross-sectional analysis of health outcomes and use of health care resources were evaluated through various metrics. Mortality rates were summarized as cases per 1000 inhabitants. In contrast, numeric health outcome variables were described by the average number of cases per person, per 100 inhabitants, or per 1000 inhabitants according to their prevalence. Average health care expenditures were reported in euro per person. Kruskal-Wallis tests, supplemented with Bonferroni-adjusted post hoc right-tailed Dunn tests, and pairwise Fisher exact tests were used to evaluate changes in the target outcomes across the risk pyramid tiers. Statistical significance was determined by considering a P value of <.05 in all analyses.

The results of the longitudinal analysis on disease prevalence and on the incidence of new disease onsets of MDD and 9 mental and somatic MDD-related chronic conditions (PR>0.80) were expressed in percentages and in per thousand (‰), respectively.

All the data analyses were conducted using R (version 4.1.1; R Foundation for Statistical Computing) [43]. The MADS algorithm was fully developed and tested in the CHSS database and transferred to the other sites using an R programming executable script.

The study is reported according to the STROBE (Strengthening the Reporting of Observational Studies in Epidemiology) [23] guidelines for observational studies.

One of the first results was the characterization of the 3 study cohorts and comparison of the sociodemographic attributes of their MADS risk groups (Table 1). All the individuals were classified into distinct risk strata based on quantiles of MADS distribution within the source population, resulting in the formation of the subsequent risk pyramid: very-low-risk tier (percentile ≤50), low-risk tier (percentile 50 to percentile 80), moderate-risk tier (percentile 80 to percentile 90), high-risk tier (percentile 90 to percentile 95), and very high–risk tier (percentile >99).

It is imperative to underscore the fundamental distinctions in the cohorts under study to comprehend the inherent sociodemographic disparities across them. Specifically, the THL and UKB cohorts predominantly consist of data derived from biobanks, specifically focusing on the middle-aged and older adult population. In contrast, the CHSS cohort represents a population-based sample encompassing the entire population spectrum.

It is worth noting that a common pattern was observed among all the cohorts in the age distribution of the citizens at risk. Although the MADS is an additive morbidity grouper, it did not monotonically increase with age. Remarkably, a notable proportion of high-risk cases were observed within the age range of 40 to 60 years, when depression typically manifests for the first time on average.

A divergence in the sex distribution across the risk strata was observable and especially noticeable in the CHSS and UKB cohorts, where the morbidity burden associated with depression and its related diseases was amplified in women (P<.001). Similarly, the disability caused by depression and its comorbidities was larger in families with fewer economic resources (P<.001). Overall, the prevalence of MDD was greater in the UKB cohort than in the other cohorts. However, upon analyzing the allocation of the population with depression in the risk pyramid, a total of 57.79% (22,238/38,479) of individuals diagnosed with MDD were categorized in the “high”- and “very high” risk tiers in the CHSS cohort, whereas the proportion of individuals diagnosed with MDD who were allocated to the tip of the risk pyramid was 40.22% (920/2287) in the THL cohort and 43.78% (23,409/53,466) in the UKB cohort.

The MADS seems to provide a novel and more comprehensive understanding of the complex nature of depression-related multimorbidity. This approach recognizes that individuals with depression often experience a range of comorbid conditions that may manifest and evolve differently over time. By capturing this dynamic aspect, the MADS offers a nuanced assessment beyond a mere checklist of discrete disorders. The novelty of the MADS approach lies in its capability to serve as the first morbidity grouper that incorporates information on disease trajectories while improving the filtering of indirect disease associations using BDMMs.

In addition to capturing disease-disease associations, the MADS endeavors to gauge their impact within the system by leveraging well-established DWs. However, despite achieving success in fulfilling the study’s objectives, it is crucial to acknowledge that this approach carries inherent limitations, as will be elaborated on in the subsequent sections of this discussion.

In this investigation, we unearthed robust correlations between the MADS risk strata and the extent of deleterious impact caused by MDD and its comorbid conditions. Such associations indicate the presence of specific health risks and an escalated use of health care resources. Furthermore, a positive association emerged between the levels of pharmacological and nonpharmacological health care expenditures and the different tiers of MADS risk. In addition, the analysis revealed an augmented risk of disease progression within the high-risk groups (high and very high risk), as indicated by a heightened incidence of new-onset depression-related illnesses within a 12-month period after the MADS assessment. Similarly, mortality rates exhibited elevated values in these high-risk groups.

The findings presented in this study are underpinned by the complementary studies conducted within the TRAJECTOME project [30] that have established a better understanding of the complex multimorbidity landscape associated with MDD across an individual’s life span, encompassing modifiable and genetic risk factors.

Despite meeting expectations and validating the hypothesis through which the study was conceived, the authors acknowledge a series of limitations leading to suboptimal results and limited potential for adaptation and generalization that should be undertaken to bring the MADS, or an indicator derived from it, to short-term real-world implementation.

In this research, the use of estimations of mean DW [44] to assess the burden of disease conditions achieved desirable results and was conceptually justified, but it undoubtedly exhibited significant limitations. In an ideal clinical scenario, each disease diagnosis indicated in the patient’s electronic medical record should be characterized by three key dimensions: (1) severity of the diagnosis, (2) rate of disease progression, and (3) impact on disability. However, the degree of maturity for characterizing the last 2 dimensions—disease progression and disability—is rather poor because of the complexities involved in their assessment. In other words, the authors acknowledge the weakness associated with the current use of DW. However, they stress the importance of incorporating such dimensions in future evolutions of the MADS.

A noteworthy aspect that should be acknowledged is that factors such as the advancements in diagnostic techniques, the digitization of medical records, and the modifications in disease taxonomy and classification over time have contributed to a more exhaustive documentation of the disease states in the most recent health records. Consequently, this fact could lead to imprecisions in estimating the disease onset ages in older individuals.

The results reported in this study not only reaffirm the well-established link between multimorbidity and adverse outcomes, such as a decline in functional status, compromised quality of life, and increased mortality rates [45], but also shed light on the significant burden imposed on individuals and health care systems. From the population-based HRA perspective, the strain on resource allocation and overall health care spending is a pressing concern that necessitates effective strategies for addressing and managing multimorbidity [46]. In this context, assessing individual health risks and patient stratification emerge as crucial approaches that enable the implementation of predictive and preventive measures in health care.

While population-based HRA tools such as Adjusted Clinical Groups, Clinical Risk Groups, or AMG have traditionally addressed this aspect, the MADS is designed to complement rather than replace those tools. This study aimed to test a method to refine existing HRA tools by aligning them with the principles of network medicine, thereby merging traditional HRA with the practical application of network medicine insights. This innovative approach holds the promise of unlocking new potential advantages and capabilities.

The strength of the MADS approach lies in using disease-disease associations drawn from the analysis of temporal occurrence patterns among concurrent diseases. This virtue allows the MADS to refine the analysis of the morbidity burden by focusing on clusters of correlated diseases, which in turn can aid in developing more tailored epidemiological risk-related studies. This refined analysis might also assist in resource allocation and inform health care policies for targeted patient groups with specific needs. Moreover, this approach holds promise for potential extrapolation to other noncommunicable disease clusters such as diabetes, cardiovascular ailments, respiratory diseases, or cancer. By leveraging this targeted approach, the MADS can be adapted to other disease clusters with shared characteristics, enabling a more precise assessment of disease burden and comorbidity patterns and thereby generating multiple disease-specific indexes.

Notably, when considering information derived from disease co-occurrence patterns, the presence or absence of certain diseases seems to correlate with the risk of developing related comorbid conditions, as elucidated in Multimedia Appendix 2. This highlights the potential for a nuanced understanding of disease relationships and their impacts on health outcomes and to implement preventive interventions to mitigate their effect. Moreover, the findings of this study highlight the potential of preventive strategies targeted at mental disorders, including substance abuse disorders, depressive disorders, and schizophrenia, to reduce the incidence of negative clinical outcomes in somatic health conditions. These important implications for clinical practice call for a comprehensive and interdisciplinary approach that bridges the gap between psychiatric and somatic medicine. By developing cross-specialty preventive strategies, health care professionals can provide more holistic and effective care for individuals with complex health needs, ensuring that their mental and physical health are adequately addressed [47].

This study provided good prospects of using disease trajectories to enhance the performance of existing state-of-the-art morbidity groupers such as AMG. Recognized for its transferability across EU regions by the EU Joint Action on implementation of digitally enabled integrated person-centered care [48], AMG stands out due to its stratification capabilities, adaptability, and distribution as open-source software, providing several advantages over its commercial counterparts. The AMG system uses disease-specific weighting derived from statistical analysis incorporating mortality and health care service use data. This method addresses the primary drawback identified in the MADS approach inherent to the use of DW while enabling the development of adaptable tools that align with the unique characteristics of each health care system. Consequently, it allows for the adjustment to the impact of specific disease conditions within distinct regions and enhances the overall applicability and adaptability of the tool. In this regard, this study offered promising insights aligned with the developers’ envisioned future features for integration into the AMG system. Serving as a proof of concept, it highlighted the potential improvements achievable within AMG by leveraging disease-disease associations, thereby shaping the road map for further AMG development.

By assessing whether the MADS is appropriate for the stratification of depression-related multimorbidity, we attempted to confirm its potential for contributing to precision medicine [49]. In the clinical arena, identifying individuals at elevated risk and customizing interventions enable health care providers to intervene proactively, potentially preventing or lessening disease progression and enhancing patient outcomes. These strategies not only yield immediate value in terms of improved patient care but also lay the foundation for the broader adoption of integrated care and precision medicine, particularly in the management of chronic conditions [50].

Incorporating systems medicine [51] methodologies and ITs has prompted significant shifts in clinical research and practice, paving the way for holistic approaches, computational modeling, and predictive tools in clinical medicine. These advancements are driving the adoption of clinical decision support systems, which use patient-specific data to generate assessments or recommendations, aiding clinicians in making informed decisions. It is well established that, to improve predictive precision and aid clinical decision-making, implementing comprehensive methodologies that consider various influencing factors from multiple sources in patient health could enhance individual prognosis estimations [52].

This integration might facilitate predictive modeling methodologies for personalized risk prediction and intervention planning. This approach, known as multisource clinical predictive modeling [53,54], enables the integration of (1) health care data and health determinants from other domains, including (2) population health registry data; (3) informal care data (including patients’ self-tracking data, lifestyles, environmental and behavioral aspects, and sensors); and, ideally, (4) biomedical research omics data. In this paradigm, it is crucial to acknowledge the pivotal role that multimorbidity groupers play in capturing the clinical complexity of individuals. Previous research [53,54] has highlighted the synergy between patient clinical complexity (eg, AMG) and acute episode severity, correlating with higher risks of adverse health events. This opens avenues for further research, exploring how adjusted morbidity indicators such as the MADS can significantly contribute to predictive modeling, aiming at supporting the implementation of cost-effective, patient-centered preventive measures to manage patients with chronic diseases and potentially delay or prevent their progression to the highest-risk levels in the stratification pyramid [55].

The MADS showed to be a promising approach to estimate multimorbidity-adjusted risk of disease progression and measure MDD’s impact on individuals and health care systems, which could be tested in other diseases. The novelty of the MADS approach lies in its unique capability to incorporate disease trajectories, providing a comprehensive understanding of depression-related morbidity burden. In this regard, the BDMM method played a crucial role in isolating and identifying true direct disease associations. The results of this study pave the way for the development of innovative digital tools to support advanced HRA strategies. Nevertheless, clinical validation is imperative before considering the widespread adoption of the MADS.