We conducted a retrospective cohort study using EMR data from both private and public data sets. Age, sex, and ICD codes were extracted. The private data set, covering 2018-2022, was obtained from a tertiary hospital in Guangxi, China, and included 21,771 hospitalizations and 6405 distinct ICD-10 (10th Revision) codes. Patients with TD were identified by ICD-10 codes E00-E07 (benign TD [BTD]) or C73 (TC). We treated each hospitalization record as an independent individual, considering that patients may have different health states and disease progression during different hospitalizations. The data cleaning process to identify the study cohort was as follows: (1) included adult patients aged ≥18 years; (2) excluded diagnostic codes from chapters XV to XXII, which pertain to pregnancy, childbirth, puerperium, perinatal conditions, and general symptoms; (3) excluded patient records with fewer than 2 ICD-10 codes; and (4) excluded rare diseases with a prevalence <2%. Ultimately, the study population comprised 18,311 patients in the private data set (Figure 1A).

Additionally, we validated our study using the Medical Information Mart for Intensive Care (MIMIC)-IV data set, a public clinical database containing over 73,000 admissions from 2008 to 2019. Diseases were initially coded using ICD-10 and ICD-9 (9th Revision) codes. To ensure consistency, we standardized all disease codes to 3-digit ICD-10 codes. Ultimately, we extracted data for 52,108 patients with TD from the MIMIC-IV data set (Figure 1B).

The original data collection for this study was approved by the Medical Ethics Committee of the First Affiliated Hospital of Guangxi Medical University (approval number 2023-E592-01), and the Ethics Committee waived the requirement for a consent form. Informed consent from patients to use their data could not be obtained because the study is retrospective, observational, and population based. The analysis of the information was anonymized, with no variables that could identify individual patients. Additionally, we collected public clinical data from the MIMIC database after completing the web course (certification number 57439457). Given that the MIMIC-IV database consists of deidentified data and is publicly accessible, and as the collection of patient information and the creation of the research resource were reviewed by the Institutional Review Board at Beth Israel Deaconess Medical Center, which granted a waiver of informed consent and approved the data-sharing initiative [14], this study was not subject to specific ethical review.

We stratified the population by sex (male or female) and disease type (TC or BTD). The prevalence of all diseases with ICD-10 codes and the 95% CI were calculated. We excluded diseases with a prevalence of <2% (P<.03, 1-sided test) in any subgroup to avoid including rare diseases. The lists of 72 and 107 diseases from the private and MIMIC-IV data sets are presented in Tables S1 and S2 in Multimedia Appendix 1, respectively. TDs, including TC (C73), nontoxic goiter (E04), hypothyroidism (E02 and E03), hyperthyroidism (E05), and other BTDs (E06), were excluded from the list of diseases, resulting in a final count of 67 and 102 comorbidities in the private and MIMIC-IV data sets, respectively. After comparing the comorbidity differences between the TC and BTD subgroups, we analyzed the comorbidity differences between patients with TC and those with specific typical BTDs (nontoxic goiter, E04; hypothyroidism, E02 and E03; and hyperthyroidism, E05). Subgroup details are presented in Figure S1 in Multimedia Appendix 1.

To evaluate the differences in comorbidity prevalence by sex and disease type, we calculated relative differences and used a Z-test for significance. The calculation for the relative difference by sex is based on equation (1) [10]:

where P_male and P_female represent the prevalences of comorbidities in males and females, respectively.

We considered comorbidity prevalence to be significantly different if the relative difference exceeded 0.1 [10] and the absolute difference was statistically significant. Among these, comorbidities with a ≥0.5-fold prevalence increase were deemed enriched. We used Spearman correlation to assess comorbidity prevalence trends with age and applied K-means clustering to group comorbidities by age-related prevalence.

We constructed a PCN to capture the coexistence of multiple diseases. In the PCN, nodes represent disease codes (ICD-10 codes at 3 digits) that are connected by edges. Node sizes are proportional to disease prevalence, and node colors represent the ICD-10 categories. The cosine index was applied to quantify the comorbidity strength of coexisting diseases, taking into account the co-occurrence and prevalence of comorbidities, thereby minimizing the influence of sample size. We defined a cutoff value to detect comorbidity coexistence measured by the cosine index by assessing the relationship between the Pearson correlation and cosine index, where the number of significant coexisting comorbidities was equal in both networks as measured by the cosine index and Pearson correlation. The cosine index and Pearson correlation coefficient are defined in equations (2) and (3), respectively, and the significance of was determined by performing a t test (2-tailed), calculated according to equation (4) [10].

where N is the total number of patients with TD; and n_a, n_b, and n_ab are the number of patients with disease a, disease b, and both diseases, respectively.

Four structural properties of the PCN were measured using network indices, including network density, degree, average degree of neighbors, and betweenness centrality [15]. Network density measures the compactness of a network by calculating the ratio of significant connections to all potential ones. A denser network indicates more connections among comorbidities. The degree of the PCN reflects its number of connections with other comorbidities. When a comorbidity is directly connected to others, those are referred to as neighbors. We calculated the average degree of neighbors to measure neighbor connectivity. Betweenness centrality represents the number of shortest paths between any 2 comorbidities. A high betweenness centrality indicates a greater likelihood of forming bridges between other comorbidities or serving as endpoints for many comorbidities [10].

To identify the most important comorbidity in the PCN, we applied the PageRank algorithm, which considers edge weights [16]. The PageRank algorithm calculates the importance of each node in a network by assigning ranks. Nodes that are connected to other nodes with a high rank receive a higher weight and are considered more central, while a higher PageRank value for a node indicates greater influence within the network. Comorbidities with the top 5 PageRank values were defined as the most important comorbidities in the PCN.

We constructed 4 separate PCNs for males, females, patients with TC, and patients with BTD. We then compared comorbidity strengths to measure disparities by sex and disease type. When a coexisting disease was unique or enriched (with a comorbidity strength of at least 0.05 higher than that of another) in a given subgroup, the difference was considered significant and was defined as an abundant connection [10].

Analyses and visualizations were performed using R 4.3.1 (R Foundation) and Python 3.7 (Python Foundation), leveraging the igraph and ggraph libraries for network visualization and property computation in R, and the pyecharts package for generating Sankey diagrams in Python.

In the private data set, the median age of patients with TD was 58 years, with 61.09% (11,187/18,311) being female and 6.25% (1144/18,311) having TC (Table 1). By contrast, in the MIMIC-IV data set, patients with TD were older, with a median age of 67 years; among them, 71.74% (37,384/52,108) were female and less than 0.93% (482/52,108) had TC. Patients with TD exhibited a high comorbidity burden, with a median number of comorbidities of 5 in the private data set and 7 in the MIMIC-IV data set. Males had more comorbidities than females, with the median number of comorbidities significantly higher in males than in females (6 vs 4, P<.001, and 8 vs 7, P<.001, in the private and MIMIC-IV data sets, respectively). Patients with TC were younger and had fewer comorbidities than patients with BTD. The median age of patients with TC was 16 and 11 years younger than that of patients with BTD, and the median number of comorbidities in patients with TC was 2 and 3 fewer than that in patients with BTD, in the private and MIMIC-IV data sets, respectively.

The mean number of comorbidities increased with age in both males and females across both data sets (Figure 2A and 2C). In the private data set, the percentage of patients with TD decreased as the number of comorbidities increased, both overall and among females. By contrast, in the MIMIC-IV data set, the percentage of patients exhibited a trend of initially increasing and then decreasing with the number of comorbidities (Figure 2B and 2D).

The incidence of patients with different types of TDs varied between the 2 data sets (Figure 3). In the private data set, nontoxic goiter was the most common BTD, affecting 67.77% (12,410/18,311) of patients, whereas in the MIMIC-IV data set, hypothyroidism was prevalent in 86.36% (44,999/52,108) of patients. In most patient subgroups, the percentage of patients decreased with the number of comorbidities in the private data set, while in the MIMIC-IV data set, the percentage of patients initially rose and then declined as the number of comorbidities increased (Figure 3B and 3E). However, across all TD subgroups within both data sets, the mean number of comorbidities increased with age, with Spearman correlations of 0.993 (P=.001) and 0.955 (P<.001) in the private and MIMIC-IV data sets, respectively (Figure 3C and 3F).

Figure 4 displays the top 15 diseases with the highest prevalence. In the private data set, the most prevalent comorbidities were hypertension, liver disease, and diabetes, with prevalences of 35.14% (6435/18,311), 21.97% (4023/18,311), and 16.82% (3080/18,311), respectively. In the MIMIC-IV data set, hypertension and diabetes were also among the most common comorbidities, with prevalence rates of 42.13% (21,951/52,108) and 20.39% (10,627/52,108), respectively. Additionally, disorders of lipoprotein metabolism (E78) showed high prevalence, at 15.29% (2800/18,311) in the private data set and 44.37% (23,120/52,108) in the MIMIC-IV data set.

Among comorbidities with significant differences (P<.05) in prevalence, the differences were substantial, with some comorbidities enriched in specific subgroups (Table 2). In the private data set, diabetes (E11), cerebrovascular diseases (I63 and I69), heart diseases (I11, I20, I25, I48, and I50), chronic obstructive pulmonary disease (J43 and J44), renal diseases (N04, N18, N20, and N28), and gout (M10) were enriched in males. Conversely, systemic connective tissue disorders (M35), dorsopathies (M48 and M50), and disorders of bone density and structure (M80 and M81) were enriched in females.

Most comorbidities had higher prevalence and were even enriched in patients with BTD, including heart disease (I20, I25, and I50), cerebrovascular disease (I63, I65, and I67), hypertension (I10), diabetes (E11), liver disease (K74), gastritis and duodenitis (K29), and renal disease (N18, N20, and N08). Notably, malignant neoplasms were enriched in the TC group, including nasopharyngeal carcinoma (C11), bronchus and lung cancer (C34), and lymph node cancers (C77). Tables S3 and S4 in Multimedia Appendix 1 present the complete prevalence differences of enriched comorbidities in the subgroups of patients with TD from the private and MIMIC-IV data sets, respectively.

Age-specific comorbidities were clustered into 3 groups in the private data set and 4 groups in the MIMIC-IV data set (Figure 5). Most comorbidities’ prevalence generally increased with age, although the rates varied by cluster. In the private data set, cluster 1 comprised 55 comorbidities with low, mostly stable prevalence. Cluster 2 included 11 comorbidities with moderate prevalence and high growth rates, such as atherosclerosis (I70), ischemic heart disease (I25), and diabetes (E11). Cluster 3 contained a single comorbidity (hypertension, I10) with the highest growth rate with age. In the MIMIC-IV data set, cluster 1 included 2 comorbidities (lipoprotein metabolism disorder, E78, and hypertension, I10) that showed significant age-related increases. Cluster 2 had 9 comorbidities with moderate increases, while clusters 3 and 4 comprised 22 and 70 comorbidities with low prevalence, respectively.

In the private data set, the PCN of patients with TD contained 72 diseases and 492 disease pairs (Figure 6A). Each disease shared a median of 13.7 significant correlations with other comorbidities. Nontoxic goiter (E04), hypothyroidism (E03), and TC (C73) were prevalent TDs, with degrees of 50, 57, and 22, respectively (Figure 6D). Hypothyroidism (E03) was strongly linked to hypertension (I10), bacterial pneumonia (J15), heart failure (I50), ischemic heart disease (I25), and anemias (D64), with cosine indexes ranging from 0.23 to 0.36 (Figure 6A and 6E). By contrast, nontoxic goiter (E04) was associated with hypertension (I10; cosine index=0.47), liver disease (K76; cosine index=0.44), atherosclerosis (I70; cosine index=0.37), gastritis and duodenitis (K29; cosine index=0.33), and disorders of the renal and ureter (N28; cosine index=0.33), differing from those of hypothyroidism (E03). Nontoxic goiter (E04) and secondary lymph node carcinoma (C77) showed the highest correlation with TC (C73), with cosine indexes of 0.32 and 0.27, respectively, while other comorbidities had a cosine index of less than 0.06 with TC.

In the private data set, hypertension (I10), atherosclerosis (I70), ischemic heart disease (I25), diabetes (E11), and liver diseases (K76) were the top 5 comorbidities with high degrees (range 32-52) and betweenness centrality (range 47.1-385.2), reflecting their frequent coexistence and mediating role in disease correlations. Among the top 103 comorbid disease pairs with a cosine index of 0.20 or higher, some connections were within the same disease system, such as 50 and 12 connections within circulatory system diseases (ICD-10 codes starting with I) and malignant neoplasms (ICD-10 codes starting with C). When the diseases mentioned above (ie, E03, E04, I10, I70, I25, E11, and K76), which account for about 10% (7/72) of the diseases in the PCN, were removed, the scale of the PCN (Figure 6C) decreased significantly, with the number of significant comorbid disease pairs dropping by 57.5% (283/492).

In the MIMIC-IV data set, there were 107 diseases with 1819 disease links in the PCN, as shown in Figure S2 in Multimedia Appendix 1. Hypertension (I10), ischemic heart disease (I25), heart failure (I50), diabetes (E10 and E11), and esophageal disease (K23) were closely connected with TDs, with the highest cosine index values of 0.60, 0.48, 0.46, 0.42, and 0.39, respectively. Figures S3 and S4 in Multimedia Appendix 1 display the PCN for each subgroup in the private and MIMIC-IV data sets, respectively.

In both data sets, the comorbidity coexistence relationships were more complex in males than in females and in patients with BTD compared with those with TC (Table 3). In the private data set, males had 60 comorbidities and females had 63, with the median numbers of connections for comorbidities of 14.0 and 11.5, respectively. Patients with TC had fewer comorbidities (lower degrees) than those with BTD. As patients aged, the PCNs exhibited an increasing average degree of nodes and their neighboring nodes across the 5 age groups, indicating increasingly complex comorbidity coexistence relationships.

Figure 7 displays the top 10 important comorbidities with the highest PageRank values in different patient groups. In both data sets, cardiovascular and cerebrovascular diseases (colored in blue, ICD-10 codes with I), such as hypertension (I10 and I11) and ischemic heart disease (I25), along with endocrine and metabolic diseases (colored in orange, ICD-10 codes with E), such as diabetes (E10 and E11) and lipoprotein metabolism disorder (E78), showed high PageRank values, indicating their significance in the PCN of patients with TD. Notably, secondary lymph node carcinoma (C77) exhibited a particularly high PageRank value within the PCN of patients with TD in the private data set, especially among patients with TC.

Figure 8 illustrates the abundant connections among patients with TD in the private data set. In the sex-specific PCN, there were 103 connections involving 55 comorbidities common to both sexes. Abundant connections in females were primarily linked to atrial fibrillation and flutter (I48), heart failure (I50), diseases of the stomach and duodenum (K29 and K31), and renal failure (N18). By contrast, abundant connections in males were mainly associated with various cardiocerebrovascular diseases (ICD-10 codes with I), malignant neoplasms (ICD-10 codes with C), bacterial pneumonia (J15), disorders of the ureter (N28), and cholelithiasis (K80). Differences in connections based on disease type were relatively minor, with 2 abundant connections in patients with TC and 5 in those without. The abundant connections in patients with TC were primarily linked to metabolic disorders (E78 and E79). By contrast, patients with BTD exhibited several abundant connections, including hypertension (I10), diabetes (E11), pneumonia (J18), and cholelithiasis (K80). The abundant connections within a specific subgroup of patients with TD from the MIMIC-IV data set are shown in Figure S5 in Multimedia Appendix 1.

In this population-based cohort study, we used a network-based approach to reveal the comprehensive comorbidity patterns in patients with TD, highlighting differences across sex, age, and TD types. Comorbidities were prevalent among patients with TD, with variations in prevalence and complex coexistence relationships influenced by sex, age, and the type of TD. Nontoxic goiter, hypothyroidism, hyperthyroidism, and TC were identified as the most common TDs. The top 5 comorbidities with the highest prevalence were hypertension, liver disease, diabetes, atherosclerosis, and disorders of lipoprotein metabolism. The number of comorbidities increased with age across all patients with TD and subgroups. Males and patients with BTD demonstrated higher rates of comorbidity, greater disease diversity, and stronger associations compared with females and patients with TC.

In this study, hypertension emerged as the most significant comorbidity, affecting 35.14% (6435/18,311) of patients with TD in the private data set. A previous study has demonstrated a positive correlation between blood pressure and thyroid-stimulating hormone levels [17]. Additionally, both blood pressure and dyslipidemia have been strongly associated with TNs [18]. Consistent with previous studies [5,19-22], we found a strong link between TD and diabetes. Specific thyroid disorders, such as TNs [22], hypothyroidism [19], thyroid dysfunction [5], and thyroid disorders [20], have been linked to diabetes. Furthermore, thyroid function significantly impacts the cardiovascular system [23]. Both hyperthyroidism and hypothyroidism are associated with increased cardiovascular morbidity and mortality, and abnormal thyroid function elevates the risk of coronary heart disease [6,24]. Our findings also revealed high PageRank values for ischemic heart diseases within the PCN of patients with TD. Additionally, thyroid hormones regulate metabolism, linking metabolic syndrome and thyroid dysfunction to significant morbidity and mortality [23]. We also observed that metabolic disorders, such as disorders of lipoprotein metabolism, were closely associated with patients with TD. Moratalla et al [8], exploring the comorbidity network of hypothyroidism, showed strong links to thyroid/larynx cancer, hyperthyroidism, anemia, and goiter. We similarly found associations with nasopharyngeal carcinoma, secondary respiratory cancers, anemia, hyperthyroidism, and nontoxic goiter. However, we revealed that hypothyroidism was strongly associated with hypertension and pneumonia, which was unreported in the study by Moratalla et al [8].

In this study, females exhibited a higher TD burden than males, especially among younger individuals. Similar findings have been widely reported in several studies [22,25]. Consistent with previous studies [22,25], TDs were more prevalent in females than males [1], and the incidence of TC was significantly higher among females, with a ratio of 2.36 in 2017 [26]. This could be attributed to factors such as physiology, pregnancy, and estrogen exposure [27,28]. Estrogen, which is higher in females, can lead to thyroid gland changes, increasing the risk of TDs [27]. Therefore, greater attention should be paid to thyroid gland examination in females.

Further, our findings revealed sex-specific comorbidity patterns in patients with TD. Males, despite lower TD incidence, tended to have more comorbidities than females, with 31 enriched disease pairs in males’ PCN compared with only 5 in females on the private data set. Moratalla et al [8] also found that comorbidity networks showed more and stronger connections in males than in females. Notably, cancers, especially malignant neoplasms of lymph nodes, bronchus and lung, and respiratory and digestive organs, demonstrated enriched connections with TC in the PCN of males. Moreover, despite evidence of higher primary TC incidence among females [26,29], males faced a greater risk of developing second primary TC compared with females [30-32]. For instance, males had a higher risk of head and neck cancer than females [30]. Our results confirmed that males were more likely to have concurrent cancers other than TC than females. Besides, males exhibited enrichment in multiple cardiocerebrovascular diseases, while females displayed fewer heart disease associations. Understanding these sex-specific comorbidity patterns aids in developing tailored treatment plans for patients with TD.

Our findings showed that older patients with TD, regardless of sex or TD type, tend to have more comorbidities, which increase with age. This aligns with previous research, including Stenholm et al [33], who reported a rise in diseases and physical difficulties with age, and Xu et al [22], who found a similar trend in TNs’ prevalence. Therefore, enhanced nursing and care for older adult patients with TD is crucial to reduce comorbidities, improve quality of life, and ease medical burdens.

Notably, the average age of patients with TC in this study was 44 years, the same for both genders. A previous study [34] also found that 45-year-old women are the most prone to papillary thyroid carcinoma, while another study [26] reported peak ages of 50-69 years for males and 15-49 years for females. However, Huang et al [2] found that the incidence of TC increased in populations aged under 40 years in China, with males showing a higher growth rate (18.6 vs 13.3). Some researchers suggested that young patients may delay seeking treatment due to underestimating the cancer risk [35]. We speculated that high workload and mental stress in young people may cause endocrine issues and hormone imbalances, contributing to TC’s younger trend. Given the link between neuropsychological functions and thyroid health [36,37], extra care is needed for young patients with TD to prevent cancer.

Besides, we observed that patients with TC, despite having a relatively low median comorbidity count (3 and 4 on the private and MIMIC-IV data sets, respectively), exhibited severe comorbidities such as cancer, heart failure, cerebral infarction, and pneumonia. This severity may stem from the aggressive nature of cancer and its tendency to metastasize, leading to life-threatening complications. Prior studies have shown higher risks of papillary TC following other neoplasms, particularly renal and breast cancers, as well as leukemias/lymphomas [38,39]. Similarly, Lian et al [40] found that subclinical hypothyroidism was common in patients with nasopharyngeal carcinoma within a year and that thyroid volume is a risk factor for radiation-induced hypothyroidism. Given these findings, heightened vigilance is necessary to prevent cancer progression and deterioration in patients with TC.

From a methodological perspective, we used cosine similarity and Pearson correlation indices to quantify node connectivity, considering disease prevalence and minimizing sample size biases. Cosine similarity effectively captured significant comorbidities in this study, as well as in previous analyses of comorbidities associated with colorectal cancer [10] and idiopathic cardiomyopathy [41]. Other similarity metrics, such as the Jaccard index [42], which quantifies the overlap between disease sets, and the observed-to-expected ratio, which considers the observed prevalence of disease pairs relative to the expected prevalence [7], were utilized in constructing cardiovascular comorbidity networks [7,42]. Additionally, various network properties were used to computationally identify key comorbidities. We applied the PageRank algorithm, which considers the number and weight of edges, the number of neighbor edges, and the centrality of neighbors, to assess the importance of comorbidities. In the comorbidity networks of both data sets in this study, the comorbidities selected by PageRank were similar, including hypertension and diabetes. PageRank has also been utilized in constructing comorbidity networks for patients with colorectal cancer [10], general patients [16], and hypothyroidism patients [8]. Additionally, other properties, such as degree [11] and degree centrality [42], have been used to identify important comorbidities and their connections.

This study has several limitations. First, apart from the MIMIC-IV data set, we collected data and conducted analyses solely in a general tertiary hospital in China, which may introduce selection bias, particularly concerning geographical and population factors. Utilizing multicenter data could yield more representative and robust results. Second, our analysis focused solely on the coexistence of comorbidities without considering their causal relationships, and many network properties remain unexplained in clinical or biological contexts. Some diseases may share common pathophysiological mechanisms or causes, leading to mutual influence. In some cases, one disease may promote the occurrence and progression of another. Therefore, further analysis is needed to explore the deeper and more detailed causal relationships between thyroid disorders and their comorbidities. Third, relying solely on ICD-10 codes for disease identification may not fully capture all cases related to thyroid disorders, potentially leading to missing information. Additionally, some diagnoses may be subject to overdiagnosis or underdiagnosis due to the complexity of clinical issues. Certain diseases may not be frequently coded, and broad categorization may obscure differences among heterogeneous conditions. Finally, we conducted only a cross-sectional analysis and did not consider the temporality and sequence of comorbidities when constructing comorbidity profiles. These patterns may change as diseases progress. By collecting more data over longer periods from multiple hospitals, a more detailed analysis could provide deeper insights into the comorbidity patterns of patients with thyroid disorders.

This data-driven exploration of comorbidity patterns among all patients with thyroid disorder and their subgroups aims to provide valuable clinical and foundational insights into the comorbidities associated with thyroid disorders. It enhances our understanding of thyroid disorders as a whole and may inform the development of more effective and integrated therapeutic strategies for patients. Additionally, the network-based methodology utilized in this study has the potential to uncover comorbidity patterns in other diseases. In the future, we plan to incorporate a wider variety of data to construct more comprehensive networks. By extracting features from these networks, we aim to leverage medical knowledge to develop machine learning–based prediction models for patients with thyroid disorders.