We considered that similar diseases are more likely to associate with similar genes (Xu and Li, 2006; Wei and Liu, 2020) and computed disease functional similarity via disease-related genes. For two diseases d_i and d_j and corresponding associated gene sets G_i = {g_i1, g_i2, …, g_ia} and G_j = {g_j1, g_j2, …, g_jb}, the functional association between gene g_k and gene set G = {g₁, g₂, …, g_l} is first defined by Eq. (2):

where FS(g_k, g_t) indicates the functional similarity between g_k and g_t by Eq. (3):

where LLS′ denotes the normalized score of LLS by Eq. (4):

where LLS represents association log-likelihood score used to evaluate the functional linkage probability between two genes provided by HumanNet (Hwang et al., 2019; Long et al., 2021), LLS_max and LLS_min denote its maximum and minimum values, respectively.

Finally, the functional similarity between d_i and d_j is computed by Eq. (5):

Microbe functional similarity matrix M_f is computed based on the method proposed by Kamneva (2017).

Based on the assumption that functionally similar diseases usually associate or disassociate with similar microbes, and disease Gaussian Interaction Profile Kernel (GIPK) similarity (Van Laarhoven et al., 2011) is computed via experimentally validated MDA network. In particular, the GIPK similarity of two diseases d_i and d_j is computed by Eq. (6):

where

and IP(d_i) denotes associations between disease d_i and each microbe, that is, the ith row of X. γ_d denotes the normalized kernel bandwidth with original bandwidth γd′ of 1, and n_d denotes the number of diseases.

Similarly, we computed the GIPK similarity matrix M_G of microbes.

We may fail to compute the functional similarity for all diseases because not all diseases have related to genes. Thus, we combined disease GIPK similarity and functional similarity by Eq. (8):

Similarly, the integrated microbe similarity S_M is computed.

For each microbe–disease pair (m_i, d_j), feature vectors of m_i and d_j are obtained based on similarity matrices S_D and S_M, respectively. Particularly, the feature vector of d_i is denoted as the similarity between d_i and all diseases. The feature vector of m_j is denoted as the similarity between m_j and all microbes. Thus, one microbe–disease pair is depicted as an (n_d+n_m)-dimensional feature vector after concatenation operation, where n_d and n_m indicate the number of diseases and microbes, respectively. In summary, there are n(n = n_d×n_m) samples (microbe–disease pairs), and each sample x_i can be represented using a d(d = n_d+n_m)-dimensional vector. For x_i, its label y_i = 1. If its corresponding microbe–disease pair is associated, otherwise y_i = 0. Consequently, an MDA matrix X with n samples is represented by Eq. (9):

The obtained features for microbe–disease pairs are highly dimensional and severely affect the classification accuracy of models. Deep learning demonstrates stronger feature learning ability than traditional dimensional reduction approaches. Thus, we designed a sparse autoencoder to reduce the feature dimensionality of each sample.

Sparse autoencoder (Andrew, 2011) is an unsupervised neural network model. It minimizes the reconstruction error and enforces sparsity constraints on all hidden nodes to obtain a more robust and meaningful representation of features and further improves the prediction performance of classification models (Makhzani and Frey, 2013). First, a high-dimensional feature vector for the microbe–disease pair is fed to an encoder by Eq. (10):

where X represents the input n samples with d-dimensional vector, H denotes the low-dimensional features after encoding, W, b, and f(·) represent the weight, bias, and encoding function of the encoder, respectively.

Next, a decoder restores the low-dimensional representation H to the same appearance as the input feature representation by Eq. (11):

where W′, b′, and g(·) represent the weight, bias, and decoding function of the decoder, respectively, and X^ denotes the learned feature representation.

To minimize the reconstruction error, we build a cost function by Eq. (12):

where λ and β denote the sparsity regularization parameter and the coefficients for L₂ regularization, respectively.

The first term MSE is mean square error. The term is used to measure the discrepancy between the input features X and the reconstructed features X^ on training data by Eq. (13):

The second term Ω_sparsity is the Kullback–Leibler divergence. The term is used to control sparsity based on the sparsity proportion ρ by Eq. (14):

where s_l and ρ^t denote the number of neurons in the lth hidden layer and the average activity of the tth neuron, respectively, KL(ρ||ρt^) denotes the relative entropy between Bernoulli random variables with mean ρ and mean ρ^t. KL(ρ||ρt^) is computed by Eq. (15):

The third term is L₂ regularization term Ω_weights. The term is used to control the weights and avoid overfitting by Eq. (16):

where n_l, s_l, and wji(l) denote the number of layers, the number of units in the lth layer, and the weight, respectively.

Each microbe–disease pair is represented as a low-dimensional vector after dimensional reduction based on a sparse autoencoder. LightGBM (Ke et al., 2017) is an optimized version of Gradient Boosting Decision Tree (GBDT) (Ye et al., 2009). It obtains better performance in the area of bioinformatics. Next, the constructed low-dimensional vector is used as the input of LightGBM (Ke et al., 2017), to classify each microbe–disease pair. For an MDA dataset D={(xi,yi)}i=1n, LightGBM intends to learn an approximation f^ to a certain function f(x) by minimizing the expectation of the loss function L(y, f(x)) by Eq. (17):

LightGBM integrates T decision trees ∑t=1Tft(X) to approximate the final model fT(X)=∑t=1Tft(X). The decision trees with J leaf nodes are expressed as w_q(x), where w_q(x) denotes the weights of all samples on leaf nodes and q(x) denotes the decision rules. Hence, The loss function of LightGBM is defined by Eq. (18):

The constant term in model (18) is removed for simplicity, and model (18) is transformed as Eq. (19):

where g_i and h_i denote the first-order and second-order derivatives of the loss function, respectively.

For a sample set, I_j related to leaf j, model (19) could be transformed as follows:

Given a tree structure q(x), the optimal leaf weight wj* of each leaf node and the maximum value of a scoring function Γ_k that evaluate the quality of q(x) are defined by Eqs. (21) and (22):

Consequently, the objective function is represented as Eq. (23):

where I_L and I_R denote the example sets on the left and right sides, respectively.

Table 1 shows accuracies, MCCs, AUCs, and AUPRs of SAELGMDA and the other four methods under CV₁. The best performance in each column is described in Tables 1–6. As shown in Table 1, SAELGMDA computed the best MCC, AUC, and AUPR on the HMDAD database and the best accuracy, MCC, AUC, and AUPR on the Disbiome database, significantly outperforming the other four MDA prediction methods under CV₁. Although accuracy was slightly less than MNNMDA and GATMDA on HMDAD, the difference was very tiny. Moreover, SAELGMDA outperformed the other methods, especially AUC and AUPR on the whole. In addition, although SAELGMDA outperformed the other four methods, all methods computed lower MCC and AUPR under CV₁, which may be caused by fewer diseases. Figure 2 shows the ROC and PR curves of the five methods on the two databases under CV₁.

Table 2 demonstrates the prediction performance of SAELGMDA and the other four methods under CV₂. The best performance in each column is described in boldface. As shown in Table 2, we observed that SAELGMDA computed the best accuracies, MCCs, and AUCs on the two databases under CV₂. In particular, SAELGMDA obtained better MCC and AUPR on the HMDAD database than ones on the Disbiome database, which may be caused by different data structures. In addition, all five MDA prediction methods computed lower MCC and AUPR on the Disbiome database. Figure 3 shows the ROC and PR curves of the five methods under CV₂.

Table 3 shows the performance of SAELGMDA and the other four methods under CV₃. The best performance in each column is described in boldface under CV₃. The results from Table 3 suggest that SAELGMDA achieved the best accuracies, MCCs, and AUCs, significantly outperforming the other four MDA prediction methods under CV₃. Moreover, the performance of all five methods under CV₃ outperforms the ones under CV₁ and CV₂, demonstrating that more samples help improve the classification performance. Figure 4 shows the ROC and PR curves of the five methods under CV₃.

To measure the MDA classification performance of LightGBM, we compared it with two classical boosting algorithms, XGBoost and NGBoost. Extreme Gradient Boosting (XGBoost) is an ensemble learning method based on a gradient boost tree and can accurately cope with multicollinearity impact and complicated non-linearity interactions (Chen and Guestrin, 2016; Zhu and Zhu, 2019). Natural Gradient Boosting (NGBoost) uses natural gradients instead of regular gradients to implement flexible probabilistic forecast (Duan et al., 2020). Tables 4–6 show the accuracy, MCC, AUC, and AUPR of LightGBM, NGBoost, and XGBoost on the Disbiome and HMDAD datasets under three cross validations. The results from Tables 4–6 indicate that LightGBM obtained better performance on the majority of conditions and can be used to improve MDA classification ability.

We compared the computational time of SAELGMDA with the other four MDA prediction models, MNNMDA, GATMDA, NTSHMDA, and LRLSHMDA. The experiments were run on a machine with an AMD EPYC 7302 CPU, a GeForce RTX 2080 Ti, and 256GB RAM on Ubuntu 20.04.4 LTS operating system. Figure 5 shows computational time (m) of the five MDA prediction models on five-fold cross validation for one time on two MDA datasets. As shown in Figure 5, SAELGMDA is the most rapid method on the HMDAD dataset and the slowest one on the Disbiome dataset. SAELGMDA need only to spend 10.57 min, although it run slowly on the Disbiome database. In summary, SAELGMDA need not too much time on the two MDA datasets.

We further predicted new MDAs based on known MDAs using SAELGMDA. The predicted top 50 MDAs are shown in Figure 6. In Figure 6, sky blue solid lines and red dotted lines represent known and unknown MDAs obtained from SAELGMDA, respectively. Deep sky blue round rectangles represent microbes and green diamonds denote diseases.

On the HMDAD database, all predicted top 50 MDAs have been known to be associated with the database. SAELGMDA predicted that Actinobacteria and liver cirrhosis have the highest association probability with the ranking of 130 among all 11,388 microbe–disease pairs. Actinobacteria have been reported to associate with liver disease (Bull-Otterson et al., 2013). The expansion of Proteobacteria and Actinobacteria has a pathogenic effect on alcoholic liver disease (Bull-Otterson et al., 2013).

In the Disbiome database, SAELGMDA predicted that Veillonella may associate with autism with a ranking of three among all 229,336 microbe–disease pairs. Zhang et al. (2018) has reported that the abundance of Veillonella was severely decreased in stools of children suffering from autism spectrum disorder. The decreasing of its abundance has been also found in subjects involved in autism (Strati et al., 2017). Furthermore, the decreased Veillonella may affect the fermentation of lactate in the autism children (Gronow et al., 2010).

Colorectal cancer is the third most frequent cause of cancer mortality worldwide, severely threatening global life and health (Biller and Schrag, 2021; Saeed et al., 2021; Wong et al., 2023). There are more than 1.85 million colorectal cancer cases and 850,000 colorectal cancer-related deaths each year. In total, 20% of patients with colorectal cancer have metastasis cancer among new colorectal cancer diagnoses. It has been reported that ~70%–75% of patients survive more than 1 year, 30%–35% more than 3 years, and fewer than 20% more than 5 years among patients diagnosed with metastatic colorectal cancer. Although colonoscopy has been widely applied to the screen, its effect on colorectal cancer remains unclear (Bretthauer et al., 2022). Table 7 shows the top 20 microbes associated with colorectal cancer on the HMDAD database.

For colorectal cancer, as shown in Table 7, 19 microbes have been confirmed to have associations with colorectal cancer by the existing literature on the top 20 inferred microbes on the HMDAD database. For example, pseudomonas is distinctly less abundant in cancer tissues than normal tissues and has been increasingly taken as an emerging clinic-related opportunistic pathogen (Decker and Palmore, 2014; Gao et al., 2015). Haemophilus parainfluenzae demonstrates higher representation in colorectal cancer subjects but is scarcely investigated in control subjects (Kasai et al., 2016). Research in 219 patients with colorectal cancer has suggested that clostridium difficile has a dense relationship with colorectal cancer (Yeom et al., 2010). Helicobacter pylori infection has been reported to be a potential risk increase factor of left-sided colorectal cancer (Zhang et al., 2012).

Moreover, we inferred that Clostridium coccoides has a possible association with colorectal cancer. Clostridium coccoides is taken as one of the most prevalent groups of bacteria in human intestines. They constitute ~60% of mucin-adhered microbiota and comprise different species with high oxygen-sensitive anaerobes (such as Clostridium, Coprococcus, Eubacterium, and Ruminococcus). They contribute to the prevention of colonization of vancomycin-resistant Enterococcus in an antibiotic-treated mouse model (Grenda et al., 2022). The association between Clostridium coccoides and colorectal cancer needs further validation.

Inflammatory bowel disease is one of the idiopathic inflammatory bowel disorders that severely affect the gastrointestinal tract. It has become a global, chronic, and life-threatening disease over the last few decades. Mak et al. (2020) predicted that patients with inflammatory bowel disease may be an exponential increase worldwide. It typically includes Crohn's disease and ulcerative colitis. It manifests progressive and unpredictable features and is partially caused by bacteria that activate patient's immune system to protect against foreign substances (Lomax et al., 2006; Kaplan and Windsor, 2021). It has a close relationship with microbes. Identification of associated microbes for the disease helps us better equip to stem its global rise in future. Table 8 lists the top 20 microbes associated with the disease on the HMDAD database.

As shown in Table 8, 19 microbes have been validated to link to inflammatory bowel disorders by existing literature on the predicted top 20 microbes associated with it on the HMDAD database. Researchers reported that Firmicutes were less represented in patients suffered from inflammatory bowel disease than healthy subjects (Sokol et al., 2009). Streptococcus and Haemophilus were highly represented in patients with inflammatory bowel disease (Heidarian et al., 2019). Prevotella was reduced in pediatric Crohn's disease (Lewis et al., 2015). Clostridium coccoides was less abundant in patients with active inflammatory bowel disease than ones in remission (Prosberg et al., 2016).

In addition, we predicted that Sphingomonadaceae dense links to inflammatory bowel disease. Sphingomonadaceae family has high abundance in marine waters, freshwater, and even drinking water. They can degrade lignin-derived compounds and refractory organic matter that comprise monocyclic and polycyclic aromatic hydrocarbons (Shen S. et al., 2022). Sphingomonadaceae are significantly accommodated to bile salts through metabolic pathways (de Vries et al., 2019). In addition, Sphingomonadaceae has a high linkage with triclosan degradation in nitrification and denitrification systems (Dai et al., 2022). Microbial communities were adapted to Bisphenol A through the selection of Sphingomonadaceae populations including Sphingobium, Novosphingobium, and Sphingopyxis. The selected Sphingomonadaceae for Bisphenol A demonstrated higher Bisphenol A metabolic activity (Oh and Choi, 2019). The association between Sphingomonadaceae and inflammatory bowel disease needs further validation.

Lung cancer is one of the leading causes of cancer-related deaths worldwide. It accounts for ~18% of global cancer deaths (Sung et al., 2021). More than 350 patients died from lung cancer each day in the United States (Siegel et al., 2022). It has the highest incidence and mortality compared with other cancer types in China (Xia et al., 2022). We used the proposed SAELGMDA model to identify potential microbes for lung cancer. Table 9 lists the top 20 microbes associated with it on the Disbiome database. As shown in Table 9, all 20 top microbes have been confirmed to be associated with lung cancer by existing literatures or the Disbiome database. The results again validated the MDA prediction performance of SAELGMDA.