The data of 2137 inpatients diagnosed with T2D were collected from the EMR database of the Second Affiliated Hospital of Shandong University of TCM from 2016 to 2021, which included all inpatient information obtained during hospitalization, such as demographic information, symptoms, laboratory or physical tests, diagnoses and treatment. Because most data were in free text that cannot be used directly for analysis, we used a clinical information extraction tool (Shu et al., 2019) to efficiently extract the biomedical entities (e.g., symptoms, diseases) from these records. Then, to normalize the various clinical term descriptions, we manually checked and standardized the terms “disease”, “botanical drug” and “drug” by referring to the 10th Revision of International Classification of Diseases (ICD-10) (World Health Organization, 2019b), the Pharmacopoeia of the People’s Republic of China 2020 Revision (ChP 2020) (Chinese Pharmacopoeia Commission, 2020), and DrugBank Online (Wishart et al., 2018), respectively. In addition, diseases with detailed ICD-10 codes were further aggregated into higher level codes. For example, the ICD-10 codes I50.903 and I50.905 were aggregated into ICD-10 code I50.9.

In this study, several external data sources were used to support this research. The efficacy of botanical drugs was extracted from ChP 2020, and human protein‒protein interactions (PPIs) were obtained from the STRING database (Szklarczyk et al., 2019). The phenotype–genotype and botanical drug–target associations were extracted from the SymMap database (Wu et al., 2019). The disease–gene associations were extracted from the MalaCards database (Rappaport et al., 2017).

Missing data and semantic mismatch were the two main challenges of EMR analysis. Therefore, we used HEMnet to address the challenges of EMR analysis by leveraging information from several external sources to supplement clinical data (Edward et al., 2017). In our study, we utilized three distinct categories of edges to create the HEMnet (Figure 1). The first two categories PPI and phenotype–genotype were drawn from the external database, while the last category was drawn directly from the EMRs.

1) PPI. This network was based on HumanNet, an external network of protein-encoding genes (Lee et al., 2011). The nodes are proteins, and the undirected edges are the interactions between proteins.

Then, HEMnet uses an embedding method, ProSNet (Wang et al., 2017), to infer relationships among its constituent nodes. ProSNet takes a heterogeneous network as input, on which it performs a novel dimensionality reduction algorithm to optimize a low-dimensional vector representation for each node. The vectors of two nodes are colocalized in the low-dimensional space if the nodes are close to each other in the heterogeneous network. After generating low-dimensional vector representations of nodes in the HEMnet, a similarity matrix was constructed according to the similarity between every two embedding vector features, which was calculated by cosine similarity. Finally, the similarity matrix was used to fill in missing features of the original patient characteristics and form the phenotypes of patients (Edward et al., 2017).

The K-means clustering (MacQueen, 1967) was used for the patient phenotype. According to the outcomes, the patients were divided into eight non-overlapping subgroups. The t-distributed stochastic neighbour embedding (t-SNE) algorithm (Cieslak et al., 2020) was used to visualize the outcomes.

The chi-square test and relative risk (RR) (Pirhaji et al., 2008; Ouimet et al., 2010) were used to assess the significance of clinical features, including symptom phenotypes, syndrome diagnoses, botanical drugs and comorbidities in eight subgroups. In this study, patients with a certain clinical feature, such as a symptom phenotype, in a particular subgroup as an exposed group, and the remaining patients with this certain clinical feature as the non-exposed group. So RR is defined as RR = (C _ij /C _i )/((C _j -C _ij )/(N-C _i )), where C _i is the number of patients in subgroup i, C _j is the number of patients with a clinical feature j, C _ij represents the number of patients in subgroup i and with a clinical feature j and N is the total number of patients in the study. A p-value <0.05, which was obtained from the chi-square test, and an RR > 1 indicated that a clinical feature was truly significant.

The GO and KEGG pathway enrichment analysis are useful to trackle the DNA-related and protein-related problems. And they offers considerable power for discovering the biological functions of genes and proteins (Chen et al., 2017). The Gene Ontology (GO) project serves as a comprehensive source for functional genomics. The project creates evidence-supported annotations to describe the biological roles of individual genome products (e.g., genes, proteins, ncRNAs, complexes) (Gene Ontology Consortium, 2015). The KEGG pathway database is the main database in Kyoto Encyclopedia of Genes and Genomes (KEGG), and it consists of manually drawn reference pathway maps together with organism-specific pathway maps (Kanehisa et al., 2017). We obtained enriched GO and KEGG pathways using the Database for Annotation, Visualization, and Integrated Discovery (DAVID), which is a web-based online bioinformatics resource that aims to provide tools for the functional interpretation of large lists of genes/proteins (Sherman et al., 2022).

As shown in the table below (Table 1), of the 2137 T2D patients, 1347 (63.03%) were male, and 790 (36.97%) were female. The ages of most T2D patients (60.60%) were between 60 and 79 years old. The average length of stay (LOS) was 14.08 ± 9.20, and for most patients (41.83%), LOS was between 8 and 14 days. We counted the distinct number of comorbidities of each patient and found that most patients had 6–10 diagnoses (56.43%).

Then, we analysed the distribution of the top five clinical features including symptom phenotypes, syndrome diagnoses, botanical drugs, and comorbidities (Table 2).

With the method introduced in the Materials and Methods, we utilized three distinct categories of edges to create the HEMnet, which contained 5,846 nodes and 125,426 connected edges. There were 3,000 symptom nodes and 2,846 gene nodes. Furthermore, there were 16,641 PPI edges, 8,749 phenotype–genotype edges, and 100,036 symptom edges.

Then, the embedding method ProSNet was used to generate low-dimensional vector representations of nodes in the HEMnet. A similarity matrix was constructed according to the similarity between every two embedding vector features, which was calculated by cosine similarity, and used to fill in missing features of the original patient characteristics to form the patient phenotypes. Finally, using the K-means clustering algorithm, eight non-overlapping patient subgroups were obtained. The t-SNE algorithm was used to visualize the clustering results (Figure 2). The numbers of patients in the eight subgroups were as follows (Table 3): H1 (n = 547, 25.60%), H2 (n = 501, 23.44%), H3 (n = 432, 20.22%), H4 (n = 298, 13.94%), H5 (n = 197, 9.22%), H6 (n = 132, 6.18%), H7 (n = 18, 0.84%), and H8 (n = 12, 0.56%).

We then selected the top 10 clinical features in these modules according to their frequency in each subgroup. Then, the RR and chi-square test (RR > 1 and p < 0.05, see Materials and methods) were used to screen the significant clinical features.

Because of fewer patients in H7 and H8 subgroups, it was less meaningful to analyse them. And since this study focused on the precision treatment of comorbidities, the H1, H2, and H4 subgroups with no significant botanical drugs and the H5 subgroup with a lower frequency of botanical drug use were excluded according to the screening results. Finally, H3 and H6 were included for further analysis.

We present the statistically significant botanical drugs, comorbidities, syndromes, and symptoms in H3 and H6 (Table 4, Table 5, Table 6, and Table 7), Poria, Astragali Radix, Glycyrrhizae Radix et Rhizoma, Cinnamomi Ramulus, and Ophiopogonis radix were the significant botanical drugs. Essential (primary) hypertension, atherosclerotic heart disease, heart failure, unstable angina, etc., were the significant comorbidities. Qi-Yin deficiency was the main significant syndrome. And chest tightness, fever, coarse lung breathing, vomiting, expectoration, etc., were the significant symptoms. In H6, Chuanxiong Rhizoma, Gastrodiae Rhizoma, and Baked Ziziphi Spinosae Semen were the significant botanical drugs. Cerebral infarction, sequelae of cerebral infarction and sequelae of intracerebral haemorrhage were the significant comorbidities. Deficient qi and blood stasis was the main significant syndrome. And poor physical activity, fever, slurring of speech, vomiting, etc., were the significant symptoms.

In this part, we explored the shared molecular associations between the significant botanical drugs and comorbidities of T2D in H3 and H6. First, we identified the distinct genes associated with each significant botanical drug and comorbidity in H3 and H6 from an external database (see Materials and methods). Then, we obtained the pathways and GO terms for the botanical drugs, comorbidities and T2D in H3 and H6 by the DAVID program (2021, see Materials and methods). Finally, we screened out pathways and GO terms with p < 0.05 from botanical drugs, comorbidities and T2D. We identified the overlapping pathways and GO terms among the botanical drugs, comorbidities, and T2D in H3 and H6 (Table 8 and Table 9). In H3, there were six overlapping pathways and eight overlapping GO terms among the botanical drugs, comorbidities, and T2D. In H6, there were no overlapping pathways among the botanical drugs, comorbidities, and T2D. Therefore, we reported on the pathways that overlapped between the two of them. There was only one overlapping GO term among the botanical drugs, comorbidities, and T2D. For example, most of the pathways and GO functions in H3 were associated with T2D, such as type II diabetes mellitus, insulin resistance, glucose metabolic process, and response to glucose. The significant botanical drugs in H3 had some overlapping pathways and GO terms with comorbidities and T2D.