Predicting microRNA-disease associations using bipartite local models and hubness-aware regression

ABSTRACT

The development and progression of numerous complex human diseases have been confirmed to be associated with microRNAs (miRNAs) by various experimental and clinical studies. Predicting potential miRNA-disease associations can help us understand the underlying molecular and cellular mechanisms of diseases and promote the development of disease treatment and diagnosis. Due to the high cost of conventional experimental verification, proposing a new computational method for miRNA-disease association prediction is an efficient and economical way. Since previous computational models ignored the hubness phenomenon, we presented a novel computational model of Bipartite Local models and Hubness-Aware Regression for MiRNA-Disease Association prediction (BLHARMDA). In this method, we first used known miRNA-disease associations to calculate the Jaccard similarity between miRNAs and between diseases, then utilized a modified kNNs model in the bipartite local model method. As a result, we effectively alleviated the detriments from ‘bad’ hubs. BLHARMDA obtained AUCs of 0.9141 and 0.8390 in the global and local leave-one-out cross validation, respectively, which outperformed most of the previous models and proved high prediction performance of BLHARMDA. Besides, the standard deviation of 0.0006 in 5-fold cross validation confirmed our model’s prediction stability and the averaged prediction accuracy of 0.9120 showed the high precision of our model. In addition, to further evaluate our model’s accuracy, we implemented BLHARMDA on three typical human diseases in three different types of case studies. As a result, 49 (Esophageal Neoplasms), 50 (Lung Neoplasms) and 50 (Carcinoma Hepatocellular) out of the top 50 related miRNAs were validated by recent experimental discoveries.

KEYWORDS:

• Microrna
• disease
• association prediction
• bipartite local models
• hubness-aware regression

Introduction

MicroRNAs (miRNAs) are one kind of short endogenous non-coding RNAs (ncRNAs) with the length of 20 ~ 25 nucleotides [1]. They can bind to the 3′ untranslated regions (UTRs) of the target messenger RNAs (mRNAs) and suppress the expression of their target mRNAs at post-transcriptional level through sequence-specific base pairing [1–4]. By this way, miRNAs can influence various biological processes including cell proliferation [5], development [6], differentiation [7], and apoptosis [8], metabolism [9,10], aging [9,10], signal transduction [11], viral infection [7] and so on. In the recent several years, thousands of miRNAs have been detected based on various experimental methods and computational models since the first two miRNAs (Caenorhabditis elegans lin-4 and let-7) were discovered more than twenty years ago [12–15]. There are 26,845 entries in the latest version of miRBase, including more than 1000 human miRNAs [16]. As accumulating experiments on miRNAs had been conducted, it was observed that miRNAs with similar sequences or secondary structures tend to play roles in similar biological processes [15]. Furthermore, the dysregulations of the miRNAs have been confirmed to be associated with the development and progression of various complex human diseases [17–19]. Recently, plenty of studies have found that miRNAs are associated with various cancers or cancer related processes [20]. For example, Chandramouli et al. [21] illustrated that there was an inverse correlation between the levels of miR-101 and EP4 receptor protein in colon cancer. The co-transfection of EP4 receptor could rescue colon cancer cells from the tumor suppressive effects of miR-101. However, EP4 receptor is negatively regulated by miR-101. Thus, the ectopic expression of miR-101 could markedly reduce the proliferation and motility of colon cancer cells. Otsubo et al. [22] demonstrated that miR-126 inhibited SRY (sex-determining region Y)-box 2 (SOX2) expression by targeting two binding sites in the 3′-UTR of SOX2 mRNA in multiple cell lines. SOX2 plays important roles in growth inhibition through cell cycle arrest and apoptosis. However, the SOX2 expression is frequently down-regulated in gastric cancers. One reason is that the highly expression of miR-126 in some cultured and primary gastric cancer cells leads to the low levels of SOX2 protein in gastric cancer cells. Besides, Nathans et al. [23] found that highly abundant miR-29a was identified in HIV-1-infected human T lymphocytes. The miR-29a specifically targets the HIV-1 3ʹ-UTR region. Specific interactions between miR-29a and HIV-1 mRNA enhance viral mRNA association with RNA-induced silencing complexes and P-body Proteins. Thus, inhibiting miR-29a enhanced HIV-1 viral production and infectivity, whereas expressing a miR-29a suppressed viral replication.

Therefore, the identifying of disease-related miRNAs could help us better understand the mechanism of complex human diseases, thus promoting the development of the diagnosis, treatment, and prevention of diseases [24–26]. However, identifying the associations between miRNAs and diseases by traditional experimental methods is demanding and costly. Therefore, finding a more economical and efficient way to predict the potential disease-related miRNAs is necessary. Today, with more and more reliable biological datasets, developing a high-efficiency computational method to uncover the potential associations between miRNAs and diseases has become a good way to overcome the drawbacks of traditional experimental methods [27–34].

In the past few years, remarkable progresses have been made in the development of prediction models for identifying potential miRNA-disease associations. Moreover, numerous computational methods have been developed based on different biological networks, systems and perspectives. Those methods could be further divided into the four categories, score function-based, machine learning-based, complex network algorithm based and multiple biological information based models [35]. Furthermore, most models were constructed based on the assumption that functionally similar miRNAs usually have connection with phenotypically similar diseases [36–38].

Lots of previous models were constructed based on network analysis algorithms. Chen et al. [30] proposed the RWRMDA model by implementing random walk with restart on the miRNA functional similarity network to make prediction for the associations. However, this model can’t make prediction for new diseases without known related miRNAs. Then, Xuan et al. [39] proposed the method named Human Disease-related MiRNA Prediction (HDMP) by integrating the known miRNA-disease associations, the miRNA functional similarity, the disease semantic similarity and the disease phenotype similarity. Furthermore, miRNAs in the same miRNA family or cluster were assigned higher weight. Then the prediction for potential miRNA-disease association was made based on the information of each miRNA’s k most similar neighbors. However, HDMP still fails to overcome the problem of making prediction for new diseases, but nonetheless it is a global method which can be used to predict all miRNA-disease associations simultaneously. Mørk et al. [40] proposed miRPD in which miRNA-Protein-Disease associations were explicitly inferred. A scoring scheme which involved the scores of miRNA-disease associations was constructed by combining the text-mining scores of miRNA-protein associations and the protein-disease associations. Xuan et al. [41] further proposed the MIRNAs associated with Diseases Prediction (MIDP) model based on random walk on the network by exploiting the more useful information from the miRNA functional similarity network, the disease semantic similarity network, and the edges between the two networks to predict reliable disease miRNA candidates. This model overcame the limitation of previous models which were unable to make association predictions for diseases without known related miRNAs by an extended walk on the network. Moreover, the negative effect of noisy data was relieved by controlling the distance the walkers go away from the labeled nodes via restarting the walking.

In order to obtain a better prediction accuracy, Chen et al. [42] proposed the Within and Between Score for MiRNA–Disease Association prediction (WBSMDA) model by integrating the known miRNA-disease associations, the miRNA functional similarity, the disease semantic similarity and the Gaussian interaction profile kernel similarity for diseases and miRNAs into a composite network. According to the network analysis, they calculated and combined the within-score and between-score to obtain the final score for potential miRNA-disease association prediction. Chen el al [43]. further presented the Heterogeneous graph inference for miRNA-disease association prediction (HGIMDA) model based on a heterogeneous graph which was constructed by the same information as WBSMDA. Then an iteration equation was constructed on the heterogeneous graph to further infer the potential miRNA-disease associations. HGIMDA could be effectively applied to new diseases and new miRNAs without any known associations. Later, Li et al. [44] proposed the Matrix Completion for MiRNA-Disease Association prediction (MCMDA) model by utilizing the matrix completion algorithm to update the adjacency matrix of known miRNA-disease associations. Through iteratively calculating the predictive scores for the candidate miRNA-disease associations, they obtained a highly reliable outcome. Then, Pasquier et al. [45] proposed the MiRAI model by representing distributional information on miRNAs and diseases in a high-dimensional vector space and defining associations between miRNAs and diseases in terms of their vector similarity in order to reveal the information attached to miRNAs and diseases by distributional semantics. Yu et al. [46] proposed the MaxFlow model based on the network information of flow model. In this method, they constructed a miRNAome-phenome network built by combining the miRNA functional similarity network, the disease semantic and phenotypic similarity network and the miRNA-disease associations network, thus uncovering the accurate miRNA-disease associations.

Additionally, there are also some computational models based on machine learning. Chen et al. [34] proposed the Regularized Least Squares for MiRNA-Disease Association (RLSMDA) method to make prediction for miRNA-disease associations. RLSMDA, a semi-supervised model, could suffice the model-training without negative samples. Then, Chen et al. [32] further proposed the Restricted Boltzmann Machine for Multiple types of MiRNA-Disease Association prediction (RBMMMDA) model which is the first model that can obtain not only new miRNA-disease associations, but also corresponding association types. This model used Restricted Boltzmann machine to predict four different types of miRNA-disease associations from a two-layered undirected miRNA-disease graph with visible and hidden units. Finally, Chen et al. [47] proposed the Ranking-based kNN for MiRNA-Disease Association prediction (RKNNMDA) model which integrated numerous information to search for k-nearest-neighbors both for miRNAs and diseases by using k-Nearest Neighbors (kNN) algorithm. In RKNNMDA, the nearest neighbors were obtained based on a descending order from the similarity scores between other miRNAs or diseases and the central miRNA or disease. Then they reranked these k-nearest-neighbors according to SVM Ranking model and implemented weighted voting to reduce the prediction bias.

The above models have their own strengths and uniqueness, but nonetheless some of them have obvious weaknesses. Moreover, although most models exhibited a sound prediction accuracy, there are still some ways to further improve those methods. To be specific, machine learning in high dimensional data spaces is particularly challenging and one of the challenges is the bad hubs. Hubness is an aspect of the curse of dimensionality. As dimensionality increases, the distribution of the number of times a point occurs among the k nearest neighbors of all other points in the dataset becomes considerably skewed to the right according to some distance measure, resulting in the emergence of hubs, that is, points which appear in many more kNN lists than other points, effectively making them ‘popular’ nearest neighbors. To be more specific, high-dimensional points that are closer to the data mean have increased the probability appearing in kNN lists of other points, even for small values of k [48]. Unfortunately, some of the hubs are bad in the sort of sense that they may mislead classification algorithms. The kNN classifier is negatively affected by the presence of ‘bad’ hubs. Although ‘bad hubs’ tend to carry more information about the location of class boundaries than other points, the ‘model’ created by the kNN classifier places the emphasis on describing non-borderline regions of the space occupied by each class. For this reason, it can be said that ‘bad’ hubs are truly bad for kNN classification, creating the need to penalize their influence on the classification decision [48]. However, the previous models failed to realize the fact that the bad hubs phenomenon which may negatively influence the prediction performance might occur under the circumstances of complex data. Furthermore, the association data between miRNAs and diseases will be widely expanded to improve the prediction accuracy in the future. Thus, there are still some challenges about how to excavate more useful information from the association data and how to overcome those undesirable problems in the complex dataset. Therefore, in this study we presented a novel model named Bipartite Local models and Hubness-Aware Regression for MiRNA-Disease Association prediction (BLHARMDA) to meet those challenges. In our model, we firstly used Jaccard-similarity to represent the similarity between the investigated miRNA (disease) and the other miRNAs (diseases) in the dataset. Then we combined it with the integrated similarity in order to describe the enhanced similarity for miRNAs and disease. In addition, it is noteworthy that combining two models to compute the semantic similarity between diseases provides a more accurate way to express the similarity between disease pairs. After that, we utilized the similarity data and miRNA-disease association data to train the Bipartite Local Model (BLM) which used an hubness-aware regression of Error Corrected k-Nearest Neighbors (ECkNN) as its local model. Consequently, we made the final prediction by utilizing this model combined with a data dimension reducing method. To evaluate our model’s effectiveness of prediction, global and local leave-one-out cross validation (LOOCV) as well as 5-fold cross validation were implemented. BLHARMDA obtained an AUC of 0.9136 and 0.8390 in global and local LOOCV respectively, which proved the high prediction accuracy of our model. Besides, the average AUC of 0.9120 ± 0.0006 in 5-fold cross validation indicated BLHARMDA’s stability in prediction. In addition, we carried out three different types of case studies on three important human diseases to further evaluate the prediction ability of BLHARMDA. We first used HMDD v2.0 database to train our model and predict for Esophageal Neoplasms. The result showed that 49 out of the top 50 potential miRNAs were validated by dbDEMC and miR2Disease [16,49]. Then we made prediction for Lung Neoplasms without any known miRNAs for this disease and 50 out of the top 50 predicted miRNAs were verified by dbDEMC, miR2Disease and HMDD v2.0. Finally, we used the known miRNA-disease associations in the HMDD v1.0 as training data and made prediction for Carcinoma Hepatocellular. As a result, 50 out of the top 50 associations were confirmed by databases and experimental literatures.

Results

Performance evaluation

We first tested the hubness in our dataset. In Figure 1, we could see that the distribution of Nk (the number of times a point occurs among the k nearest neighbors of all other points) skewed to right. In other words, few points appeared in many more kNN lists than other points. Then, to evaluate the prediction accuracy of BLHARMDA, we took advantage of the HMDD v2.0 database which offered 5430 miRNA-disease associations between 383 diseases and 495 miRNAs. We implemented global LOOCV, local LOOCV and 5-fold cross validation methods based on the known miRNA-disease associations in the HMDD v2.0. In LOOCV, we picked out one known association without repetition as the test sample at each time and considered the other known associations as training examples until all known associations were evaluated. The difference between global LOOCV and local LOOCV was that all unknown miRNA-disease pairs were considered as the candidate samples in the global LOOCV while only those miRNA-disease pairs in which the miRNAs without any known associations with the investigated disease were considered as the candidate samples in the local LOOCV. As for 5-fold cross validation, we randomly split our data set, the known miRNA-disease associations, into 5 disjoint subsets with equal size. Then we took out one subset without repetition as the test sample at each time and the rest four subsets were regarded as training samples. Besides, in 5-fold cross validation, all unknown miRNA-disease pairs were considered as candidate samples in the same way as global LOOCV. Subsequently, we implemented BLHARMDA to obtain a predicted association score matrix, and ranked the score of each test sample against the score of the candidate samples. The whole process was repeated 100 times to obtain an evaluation of BLHARMDA’s prediction accuracy. As a result, for both LOOCV and 5-fold cross validation, we obtained the ranking for each sample. Based on the rankings obtained from cross validation, the model we proposed would be considered to successfully predict an association if the ranking of a test sample was above a given threshold. Then we drew the receiver operating characteristics (ROC) curves through plotting the true positive rate (TRR sensitivity) versus the false positive rate (FRR 1-specificity) at different thresholds, and calculated the area under the ROC curves (AUC) which was an evaluation metric widely used in the description of the prediction accuracy of a computational model. Specifically, sensitivity refers to the percentage of the true positive samples whose rankings are higher than the given threshold in the whole positive samples. Meanwhile, specificity denotes the percentage of negative samples with rankings lower than the threshold in the whole negative samples. The AUC of 1 indicated that all test samples were correctly predicted, and the AUC of 0.5 indicated that the model randomly predicted the test samples. The Figure 2 illustrated that, BLHARMDA achieved an AUC of 0.9141 in Global LOOCV whose performance was superior to RLSMDA (0.8426), WBSMDA (0.8030), MCMDA (0.8749), MaxFlow (0.8624) and RKNNMDA (0.7159). In Local LOOCV, BLHARMDA achieved an AUC of 0.8390 while RLSMDA, WBSMDA, MCMDA, MaxFlow, RKNNMDA, RWRMDA, MiRAI and MIDP obtained AUCs of 0.6953, 0.8031, 0.7718, 0.7774, 0.8221, 0.7891, 0.6299 and 0.8196, respectively. RWRMDA, MiRAI and MIDP were not applicable to global LOOCV since they were based on the local method which could only be used to make prediction for one disease at a time. In our evaluation, MiRAI did not achieve the performance as in literature [45] since this method was based on collaborative filtering which did not perform well in sparse data set. Our data set was sparse in which the average associations for one disease were only 14 while the MiRAI was tested on the dataset in which one disease was associated with at least 20 miRNAs in literature [45]. Thus, we could consider that our model achieved more reliability in the prediction of miRNA-disease associations.

Figure 1. The distribution of the number of kNN lists a point showed for the miRNAs and diseases in our dataset. Nk means the number of kNN lists a point showed, P(Nk) means the number of points showed in Nk kNN lists.

Figure 2. Performance comparison between BLHARMDA and seven classical disease-miRNA association prediction models (MCMDA, HGIMDA, WBSMDA, HDMP, RLSMDA, MaxFlow and RKNNMDA) in terms of ROC curves and AUCs based on global LOOCV and comparison of AUCs based on the local LOOCV between BLHARMDA and above seven models and three local models (RLSMDA, MiRAI, MIDP). As a result, BLHARMDA outperformed other models by achieving an AUC of 0.9141 in global LOOCV and an AUC of 0.8390 in local LOOCV.

In the case of 5-fold cross evaluation, BLHARMDA obtained an averaged AUC of 0.9120+/-0.0006 which was superior to RLSMDA (0.8569+/-0.0020), WBSMDA (0.8185+/-0.0009), MCMDA (0.8767+/-0.0011), MaxFlow (0.8579+/-0.0010) and RKNNMDA (0.6723+/-0.0027). RWRMDA, MiRAI and MIDP were still not included in this evaluation since this evaluation was only applicable to global models. It was noticeable that both the average AUC and the standard deviation of AUC of BLHARMDA performed better than other models, which indicated that BLHARMDA had good stability and low predication variance.

Since the performance of our model in local LOOCV is marginally higher than those compared methods, we did a Kolmogorov-Smirnov test between our algorithm and the other ten algorithms compared in local LOOCV (see Table 1). The p-value in all test are quite low and lower than 0.01 which means there are significant difference between the results of BLHARMDA and other ten models in local LOOCV.

Table 1. Kolmogorov-Smirnov test based on the local LOOCV results between BLHARMDA and other ten algorithms compared in local LOOCV.

Algorithms	HDMP	HGIMDA	MaxFlow	MCMDA	MIDP
p-value	1.0117e-82	2.3893e-08	3.5906e-87	5.0042e-120	5.1596e-11
Algorithms	MiRAI	RKNNMDA	RLSMDA	RWRMDA	WBSMDA
p-value	2.2998e-317	2.7878e-17	6.9061e-173	2.1338e-28	1.5126e-29

Case studies

In order to evaluate the predictive ability of BLHARMDA and verify the correctness and reasonableness of the predicting outcome obtained by our model in realistic cases, we implemented 3 different types of case studies on three important human complex diseases.

In the first case study, we used the known miRNA-diseases associations in the HMDD v2.0 database as the training set of our model, then we ranked the candidate miRNAs for each disease respectively according to their predicted scores. In order to further promote experimental validation, we produced a complete prediction list for all the 383 diseases in HMDD v2.0 predicted by BLHARMDA (See Supplementary Table 1). Consequently, we verified the top 50 candidate miRNAs in two databases, dbDEMC and miR2Disease. After comparing the entries between HMDD v2.0 and dbDEMC/miR2Disease, we found that 546 of the 5430 known associations in HMDD v2.0 also existed in dbDEMC and 232 associations in HMDD v2.0 also existed in miR2Disease. However, there was no overlap between the training samples and the predicted results because in case studies only candidate miRNAs which have no known associations with the investigated disease according to HMDD v2.0 were ranked and confirmed. Thus, the detriment from the overlap could be eliminated and it avoided the situation that the miRNAs associated with the investigated disease were easier to be evaluated if there existed overlap between the training database and the validation databases. Therefore, it guaranteed that a believable validation could be obtained when we evaluated our outcome by dbDEMC/miR2Disease. Esophageal Neoplasms (EN) was investigated in our first case study. Esophageal cancer is the eighth most common incident cancer in the world and sixth in cancer mortality [50]. In the United States, 4–10 in 100,000 persons succumb to the disease per year [51]. We used BLHARMDA to predict the potential related miRNAs for EN and 49 out of the top 50 (see Table 2) and 9 out of the top 10 potential miRNAs were confirmed by dbDEMC or miR2Disease. For example, miR-21 (1st in the prediction list) targets PDCD4 at the posttranscriptional level and regulates cell proliferation and invasion in esophageal squamous cell carcinoma (ESCC) [52]; miR-17 (2nd in the prediction list) can serve as potential prognostic biomarkers in ESCC which are associated with some clinicopathologic factors [53]; Downregulated expressions of miR-155 (3rd in the prediction list) in plasma were significantly associated with increasing risks of esophageal cancer [54]. Further, since miR-155 displayed significantly lower expression in plasma of patient with esophageal cancer and serum miR-155 is associated with lymphocyte-mediated immune responses, the expression of circulating miR-155 may reflect compromised immunoreactivity in patients with esophageal cancer [55].

Table 2. Prediction of the top 50 potential Esophageal Neoplasms-related miRNAs based on known associations in HMDD v2.0 database. The first column records top 1–25 related miRNAs. The third column records the top 26–50 related miRNAs. The evidences for the associations were either dbDEMC and miR2Disease.

miRNA	Evidence	miRNA	Evidence
hsa-mir-21	dbDEMC;miR2Disease	hsa-mir-222	dbDEMC
hsa-mir-17	dbDEMC	hsa-mir-181b	dbDEMC
hsa-mir-155	dbDEMC	hsa-mir-200b	dbDEMC
hsa-mir-20a	dbDEMC	hsa-mir-31	dbDEMC
hsa-mir-146a	dbDEMC	hsa-mir-15a	dbDEMC
hsa-mir-145	dbDEMC	hsa-let-7c	dbDEMC
hsa-mir-34a	dbDEMC	hsa-mir-29c	dbDEMC
hsa-mir-125b	dbDEMC	hsa-mir-146b	dbDEMC
hsa-mir-92a	unconfirmed	hsa-mir-34c	dbDEMC
hsa-mir-126	dbDEMC	hsa-let-7e	dbDEMC
hsa-mir-221	dbDEMC	hsa-mir-200a	dbDEMC
hsa-mir-18a	dbDEMC	hsa-mir-142	dbDEMC
hsa-let-7a	dbDEMC	hsa-mir-30a	dbDEMC
hsa-mir-16	dbDEMC	hsa-mir-9	dbDEMC
hsa-mir-19b	dbDEMC	hsa-let-7d	dbDEMC
hsa-mir-143	dbDEMC	hsa-mir-182	dbDEMC
hsa-mir-200c	dbDEMC	hsa-mir-199a	dbDEMC
hsa-mir-1	dbDEMC	hsa-mir-203	dbDEMC;miR2Disease
hsa-mir-223	dbDEMC;miR2Disease	hsa-mir-7	dbDEMC
hsa-mir-210	dbDEMC	hsa-mir-106b	dbDEMC
hsa-mir-19a	dbDEMC	hsa-mir-10b	dbDEMC
hsa-let-7b	dbDEMC	hsa-mir-24	dbDEMC
hsa-mir-29a	dbDEMC	hsa-mir-205	dbDEMC;miR2Disease
hsa-mir-181a	dbDEMC	hsa-let-7g	dbDEMC
hsa-mir-29b	dbDEMC	hsa-mir-150	dbDEMC

In our second case study, we hid all known related miRNAs of each investigated disease. Then we trained the model and ranked all candidate miRNAs for each investigated disease. After that, we especially investigated the Lung Neoplasms (LN) in this case study. Lung cancer is the most common cause of cancer deaths worldwide in both man and woman [56]. According to the American Cancer Society, LN account for about 13% of all new cancers and 27% of all cancer deaths [56]. There are estimated 1.4 million deaths of lung cancer each year [56–59]. Especially, lung cancer has become the first cause of death among people with malignant tumors in China and the registered lung cancer mortality rate in China has increased by 464.84% in the past three decades [60]. The five-year survival rate of lung cancer is much lower than many other leading cancers [56,59,61–63]. We used BLHARMDA to uncover the potential related miRNAs for LN. As a result, 50 out of the top 50 (see Table 3) and 10 out of the top 10 potential related miRNAs were confirmed by dbDEMC or miR2Disease or HMDD v2.0. For example, in a study performed on 48 pairs of non-small cell lung cancer (NSCLC) specimens, the overexpression of miR-21 (1st in the prediction list) was inversely correlated with overall survival of the patients, suggesting that a high level of miR-21 is an independent negative prognostic factor for survival in NSCLC patients [64]. High expression of miR-155 (2nd in the prediction list) was correlated with poor survival of lung cancer by univariate analysis as well as multivariate analysis for miR-155. The miRNA expression signature on the outcome of a study indicated that miRNA expression profiles were diagnostic and prognostic markers of lung cancer [65]. MiR-17-5p (3rd in the prediction list) have been found to be continuously expressed in small-cell lung cancer (SCLC) cells. The overexpression of miR-17-5p may serve as a fine-tuning influence to counterbalance the generation of DNA damage in RB-inactivated SCLC cells, thus reducing excessive DNA damage to a tolerable level and consequently leading to genetic instability [66].

Table 3. Prediction of the top 50 potential Lung Neoplasms-related miRNAs based on known associations in HMDD v2.0 database. In this case study, we hided all known related miRNAs for each investigated disease before the prediction process. The first column records top 1–25 related miRNAs. The third column records the top 26–50 related miRNAs. The evidences for the associations were either HMDD v2.0, dbDEMC and miR2Disease.

miRNA	Evidence	miRNA	Evidence
hsa-mir-21	dbDEMC;miR2Disease;HMDD2	hsa-mir-31	dbDEMC;miR2Disease;HMDD2
hsa-mir-155	dbDEMC;miR2Disease;HMDD2	hsa-mir-181b	dbDEMC;HMDD2
hsa-mir-17	miR2Disease;HMDD2	hsa-mir-200b	dbDEMC;miR2Disease;HMDD2
hsa-mir-146a	dbDEMC;miR2Disease;HMDD2	hsa-mir-222	dbDEMC;HMDD2
hsa-mir-20a	dbDEMC;miR2Disease;HMDD2	hsa-mir-15a	dbDEMC
hsa-mir-145	dbDEMC;miR2Disease;HMDD2	hsa-mir-146b	miR2Disease;HMDD2
hsa-mir-34a	dbDEMC;HMDD2	hsa-let-7c	dbDEMC;miR2Disease;HMDD2
hsa-mir-125b	miR2Disease;HMDD2	hsa-mir-29c	dbDEMC;miR2Disease;HMDD2
hsa-mir-126	dbDEMC;miR2Disease;HMDD2	hsa-let-7e	miR2Disease;HMDD2
hsa-mir-92a	HMDD2	hsa-mir-142	HMDD2
hsa-mir-221	dbDEMC;HMDD2	hsa-mir-30a	miR2Disease;HMDD2
hsa-mir-18a	dbDEMC;miR2Disease;HMDD2	hsa-mir-34c	dbDEMC;HMDD2
hsa-let-7a	dbDEMC;miR2Disease;HMDD2	hsa-mir-200a	dbDEMC;miR2Disease;HMDD2
hsa-mir-16	dbDEMC;miR2Disease	hsa-mir-199a	dbDEMC;miR2Disease;HMDD2
hsa-mir-19b	dbDEMC;HMDD2	hsa-mir-9	miR2Disease;HMDD2
hsa-mir-143	dbDEMC;miR2Disease;HMDD2	hsa-let-7d	dbDEMC;miR2Disease;HMDD2
hsa-mir-1	dbDEMC;miR2Disease;HMDD2	hsa-mir-7	miR2Disease;HMDD2
hsa-mir-200c	dbDEMC;miR2Disease;HMDD2	hsa-mir-106b	dbDEMC
hsa-mir-29a	dbDEMC;miR2Disease;HMDD2	hsa-mir-203	dbDEMC;miR2Disease;HMDD2
hsa-mir-223	HMDD2	hsa-mir-182	dbDEMC;miR2Disease;HMDD2
hsa-mir-29b	dbDEMC;miR2Disease;HMDD2	hsa-mir-24	miR2Disease;HMDD2
hsa-mir-210	dbDEMC;miR2Disease;HMDD2	hsa-mir-150	dbDEMC;miR2Disease;HMDD2
hsa-let-7b	miR2Disease;HMDD2	hsa-mir-10b	dbDEMC;HMDD2
hsa-mir-19a	dbDEMC;miR2Disease;HMDD2	hsa-let-7g	dbDEMC;miR2Disease;HMDD2
hsa-mir-181a	dbDEMC;HMDD2	hsa-mir-205	dbDEMC;miR2Disease;HMDD2

Lastly, in the third case study, we trained our model based on the known miRNA-disease associations in the HMDD v1.0 and then used the model to predict the scores of the miRNAs related with each investigated disease to examine the applicability of BLHARMDA to different datasets other than in HMDD v2.0. After that we evaluated the top 50 related miRNAs for each investigated disease by HMDD v2.0, dbDEMC and miR2Disease. In particular, we investigated Hepatocellular carcinoma (HCC) in this case study. HCC is one of the most common malignancies worldwide, and is the fourth most common cause of mortality [67]. In addition, its incidence is increasing in many countries [68,69]. HCC is difficult to manage, as compared with other common malignant diseases due to the high percentage of co-existing liver cirrhosis. The impaired liver function caused by liver cirrhosis often restricts treatment options, even for early HCC. We implemented BLHARMDA to find the potential related miRNAs for HCC. As a result, all out of the top 50 (see Table 4) potential miRNAs were confirmed by dbDEMC or miR2Disease or HMDD. Among the top 10 potential miRNAs, most of them were confirmed by experiments. For example, a study [70] showed that the miR-21 (1st in the prediction list), the miR-17–92 polycistron which mainly includes miR-17-5p (2nd in the prediction list) and miR-20a (3rd in the prediction list) exhibited increased expression in HCC cell lines than those observed in normal liver.

Table 4. Prediction of the top 50 potential Carcinoma Neoplasms-related miRNAs based on known associations in the HMDD v1.0 database. The first column records top 1–25 related miRNAs. The third column records the top 26–50 related miRNAs. The evidences for the associations were either HMDD v2.0, dbDEMC and miR2Disease.

miRNA	Evidence	miRNA	Evidence
hsa-mir-21	miR2Disease;HMDD	hsa-let-7i	dbDEMC;HMDD
hsa-mir-17	miR2Disease;HMDD	hsa-mir-126	dbDEMC;miR2Disease;HMDD
hsa-mir-20a	dbDEMC;miR2Disease;HMDD	hsa-mir-29b	dbDEMC;HMDD
hsa-mir-155	dbDEMC;miR2Disease;HMDD	hsa-mir-106b	dbDEMC;miR2Disease;HMDD
hsa-mir-146a	dbDEMC;miR2Disease;HMDD	hsa-mir-143	dbDEMC;miR2Disease
hsa-mir-18a	dbDEMC;miR2Disease;HMDD	hsa-let-7f	miR2Disease;HMDD
hsa-mir-19b	miR2Disease;HMDD	hsa-let-7g	miR2Disease;HMDD
hsa-let-7a	dbDEMC;miR2Disease;HMDD	hsa-mir-181b	dbDEMC;miR2Disease;HMDD
hsa-mir-221	dbDEMC;miR2Disease;HMDD	hsa-mir-146b	HMDD
hsa-mir-19a	dbDEMC;miR2Disease;HMDD	hsa-mir-29a	dbDEMC;HMDD
hsa-mir-16	dbDEMC;miR2Disease;HMDD	hsa-mir-214	dbDEMC;miR2Disease;HMDD
hsa-mir-1	miR2Disease;HMDD	hsa-mir-141	miR2Disease;HMDD
hsa-mir-222	dbDEMC;miR2Disease;HMDD	hsa-mir-127	miR2Disease;HMDD
hsa-mir-92a	miR2Disease;HMDD	hsa-mir-9	miR2Disease
hsa-let-7e	dbDEMC;miR2Disease;HMDD	hsa-mir-132	miR2Disease
hsa-mir-145	dbDEMC;miR2Disease;HMDD	hsa-mir-200a	dbDEMC;miR2Disease;HMDD
hsa-mir-223	miR2Disease;HMDD	hsa-mir-106a	dbDEMC;miR2Disease;HMDD
hsa-let-7b	miR2Disease;HMDD	hsa-mir-133a	miR2Disease
hsa-mir-15a	dbDEMC;miR2Disease;HMDD	hsa-mir-29c	dbDEMC;HMDD
hsa-mir-125b	miR2Disease;HMDD	hsa-mir-150	dbDEMC;miR2Disease;HMDD
hsa-mir-200b	miR2Disease;HMDD	hsa-mir-125a	dbDEMC;miR2Disease;HMDD
hsa-let-7d	miR2Disease;HMDD	hsa-mir-24	miR2Disease;HMDD
hsa-mir-34a	dbDEMC;miR2Disease;HMDD	hsa-mir-30d	dbDEMC;HMDD
hsa-let-7c	dbDEMC;miR2Disease;HMDD	hsa-mir-34c	HMDD
hsa-mir-199a	dbDEMC;miR2Disease;HMDD	hsa-mir-20b	dbDEMC;HMDD

Moreover, since the HMDD v2.0 database used in our research was released in August 2013, we added a Table for the case study 1 by using the studies published after 2014 in PubMed to further evaluate the top 50 potential Esophageal Neoplasms-related miRNAs. Here, we only collected the literatures after 2014 in order to provide a fairer evaluation. The result showed that 46 out of the top 50 predictions were confirmed by the experimental literatures published after 2014 in PubMed (see Table 5).

Table 5. Here, we used the studies published after 2014 in PubMed to further evaluate the top 50 potential Esophageal Neoplasms-related miRNAs. We provided the PMID and published year of these studies in the table. The result showed that 46 out of the top 50 predictions were confirmed by the experimental literatures published after 2014 in PubMed.

miRNA	Evidence (PMID)	Publication time	miRNA	Evidence (PMID)	Publication time
hsa-mir-21	29,568,234	2018 Mar	hsa-mir-222	unconfirmed
hsa-mir-17	28,002,789	2017 Feb	hsa-mir-181b	27,189,061	2016 May
hsa-mir-155	29,660,336	2018 Jun	hsa-mir-200b	27,496,804	2016 Sep
hsa-mir-20a	27,508,097	2016 Jul	hsa-mir-31	25,568,668	2014 Nov
hsa-mir-146a	27,832,663	2016 Dec	hsa-mir-15a	27,802,201	2016 Sep
hsa-mir-145	29,852,786	2018 May	hsa-let-7c	unconfirmed
hsa-mir-34a	29,094,237	2018 Jun	hsa-mir-29c	25,928,282	2015 Apr
hsa-mir-125b	29,749,531	2018 Jul	hsa-mir-146b	24,589,738	2014 Mar
hsa-mir-92a	25,826,212	2015 Mar	hsa-mir-34c	28,852,310	2017 Aug
hsa-mir-126	28,536,606	2016 Jul	hsa-let-7e	28,408,353	2017 Jul
hsa-mir-221	27,501,171	2016 Nov	hsa-mir-200a	28,025,999	2016 Nov
hsa-mir-18a	27,291,152	2016 Aug	hsa-mir-142	unconfirmed
hsa-let-7a	29,393,461	2018 Apr	hsa-mir-30a	29,259,372	2017 Dec
hsa-mir-16	24,852,767	2014 May	hsa-mir-9	25,375,090	2014 Nov
hsa-mir-19b	25,117,812	2014 Oct	hsa-let-7d	unconfirmed
hsa-mir-143	26,427,659	2016 Mar	hsa-mir-182	26,498,375	2015 Oct
hsa-mir-200c	29,113,666	2017 Dec	hsa-mir-199a	26,717,044	2016 Feb
hsa-mir-1	26,414,725	2016 Feb	hsa-mir-203	25,216,463	2015 Apr
hsa-mir-223	24,390,317	2014 Mar	hsa-mir-7	29,906,417	2018 Jun
hsa-mir-210	28,968,550	2017 Nov	hsa-mir-106b	27,619,676	2016 Nov
hsa-mir-19a	28,621,611	2017 May	hsa-mir-10b	26,554,762	2015 Nov
hsa-let-7b	24,576,011	2014 May	hsa-mir-24	25,591,590	2015 Jan
hsa-mir-29a	25,435,940	2015 Jan	hsa-mir-205	27,974,696	2017 Jan
hsa-mir-181a	25,230,784	2014 Nov	hsa-let-7g	26,655,271	2016 Feb
hsa-mir-29b	25,866,219	2015 Apr	hsa-mir-150	29,081,413	2017 Dec

Discussion

The investigating of the miRNA-disease associations has produced plenty of potential diagnostic methods since an increasing number of models have been proposed to uncover the relationship between miRNAs and diseases. In this paper, we proposed a novel computational model based on bipartite local models with error corrected k-nearest neighbors as the local model to predict the potential miRNA-disease associations. BLHARMDA firstly processed the data by integrating the miRNA functional similarity or disease semantic similarity with the miRNA or disease Gaussian interaction profile kernel similarity to obtain the integrated miRNA or disease similarity. Then we calculated the Jaccard-similarity for miRNAs and diseases and combined the integrated similarity and the Jaccard-similarity to construct the enhanced similarity-based representation. After that, we used the BLM to predict the associations to obtain our final prediction. Moreover, we carried out cross validation for BLHARMDA, and obtained an AUC of 0.9141 in Global LOOCV which outperformed some previous models such as RLSMDA, WBSMDA, MCMDA, MaxFlow and RKNNMDA. In Local LOOCV, BLHARMDA obtained an AUC of 0.8390 which still outperforms the above models as well as the RWRMDA, MiRAI and MIDP. Then we used the 5-fold cross validation to evaluate the prediction stability of BLHARMDA. The low standard derivation as well as the high mean value in the result showed high stability and the high accuracy of BLHARMDA. Furthermore, we carried out three types of case study to demonstrate the prediction accuracy of BLHARMDA. As a result, a majority of the top 50 potential related miRNAs for each disease were confirmed by the databases or experimental literatures.

The high precision and stability of BLHARMDA are based on three factors. Firstly, we made use of both the intrinsic similarity of miRNAs and diseases and the association information between miRNAs and diseases by combining the integrated similarity and the Jaccard similarity. Secondly, we used an hubness-aware regression as the local model in BLM which enables the model to effectively alleviate the negative influence from the presence of the bad-hubs. Since bad hubs are expected in the complex data such as the miRNA-disease association database we used in this study, overcoming the bad-hubs problem can help us improve the prediction accuracy of the model. Last, we used the experimentally confirmed miRNA-disease associations from the highly reliable HMDD v2.0 database to train our model to predict potential associations between miRNAs and diseases. Furthermore, plenty of biological information database such as the functional similarity between miRNAs and the Directed Acyclic Graph of the diseases were utilized to improve the prediction accuracy. By combining biological knowledge, the model could obtain a better accuracy since the biological information could help us make use of some intrinsic differences and connections involved in dataset which are helpful in solving some biological problems. Moreover, nowadays biologists have produced massive research results and brought about a lot of biological information datasets which are convenient to be used in the models.

However, there are still some limitations in BLHARMDA which need to be improved in the future. Firstly, the miRNA and disease enhanced similarity representation in this study may not be the perfect similarity calculation method. The Jaccard similarity we used to represent the interaction information between miRNAs and diseases can still excavate the connection between them at a basic level. Developing a better similarity calculation method to make use of the association information in a more efficient way may produce a better outcome from this model. Secondly, the local model in BLM in this study is based on k-nearest neighbors regression. The performance of the algorithm might be further improved by using other hubness-aware regression models based on support vector machines, neural network, etc. Thirdly, integrating more information about the miRNAs and diseases to the input data could also facilitate the prediction since we only used the miRNA’s functional similarity and disease’s semantic similarity in this study. Therefore, BLHARMDA could still be developed better in the future.

Materials and methods

Human miRNA-disease associations

Recently, more and more miRNA-disease associations have been discovered by biological experiments. In this paper, the known miRNA-disease associations dataset was acquired from HMDD v2.0 which contained 5430 associations between 495 miRNAs and 383 diseases [71]. We further constructed an adjacency matrix A to store the information of known miRNAs-disease associations from HMDD v2.0. In the adjacency matrix, A(i, j) equal to 1 means the i-th miRNA $m_i$$m_i$ is related to the j-th disease$d_j$$d_j$, otherwise A(i, j) equal to 0. Furthermore, we used nm to denote the number of miRNAs and nd denote the number of diseases.

MiRNA functional similarity

We attained the miRNA functional similarity scores from http://www.cuilab.cn/files/images/cuilab/misim.zip. Wang et al. [72] proposed the calculation method of miRNA functional similarity based on the assumption that the miRNAs with a high functional similarity are more likely to correlate with diseases with a high phenotypical similarity. Thus, we could construct the $nm\times nm$$nm\times nm$ functional similarity matrix FS. The FS(i, j) entry ranging from zero to one denotes the functional similarity score between $m_i$$m_i$ and$m_j$$m_j$.

Disease MeSH descriptors and directed acyclic graph

According to the U.S. National Library of Medicine (MeSH) at https://www.nlm.nih.gov/mesh/[39], we constructed a Directed Acyclic Graph (DAG) to describe the semantic information of a disease$d_i$$d_i$. The DAG obtained from MeSH provided a strict system for disease classification for the research of the relationship among various diseases [73]. MeSH descriptors included 16 categories: Category A for anatomic terms, Category B for organisms, Category C for diseases, Category D for drugs and chemicals and so on. We used MeSH descriptor of Category C for each disease in this paper.

The nodes in DAG represent disease MeSH descriptors and all the MeSH descriptors in the DAG are connected by a direct edge from a more general term (parent node) to a more specific term (child node). Each MeSH descriptor has one or more tree numbers to numerically define its location in the DAG. The tree numbers of a child node are defined as the codes of its parent nodes appended by the child’s information. For the disease$d_i$$d_i$, its DAG was defined as DAG($d_i$$d_i$) = ($d_i,$$d_i,$D($d_i$$d_i$),E($d_i$$d_i$)), where D($d_i$$d_i$) denotes the node set include the $d_i$$d_i$ and its ancestor diseases and E($d_i$$d_i$) is the set of corresponding direct edges from a parent node to a child node which represent the relationship between these two nodes.

Then, we used DAG to construct the disease semantic similarity matrix. Since there are two models to calculate the semantic similarity between disease pairs, we used both the two models and average them as the semantic similarity matrix in this study.

Disease semantic similarity 1

In the first disease semantic similarity model, as described in the former algorithm [72], the contribution of disease t in $DAG(d_i)$$DAG(d_i)$ to the semantic value of $d_i$$d_i$ and the semantic value of $d_i$$d_i$ could be computed by the following equations respectively.

(1) $D{1}_{d_i}\left(t\right)=\left\{\begin{array}{c}1ift=d_i\\ max\left\{\mathrm{\Delta }\ast D{1}_{d_i}\left(d^{\prime }\right)|d^{\prime }\in childrenoft\right\}ift\ne d_i\end{array}\right.$$D{1}_{d_i}\left(t\right)=\left\{\begin{array}{c}1ift=d_i\\ max\left\{\mathrm{\Delta }\ast D{1}_{d_i}\left(d^{\prime }\right)|d^{\prime }\in childrenoft\right\}ift\ne d_i\end{array}\right.$(1)

(2) $DV1\left(d_i\right)=\sum _{d\in D\left(d_i\right)}D{1}_{d_i}\left(d\right)$$DV1\left(d_i\right)=\sum _{d\in D\left(d_i\right)}D{1}_{d_i}\left(d\right)$(2)

Here, $\mathrm{\Delta }$$\mathrm{\Delta }$is the semantic contribution delay factor, and this factor was set to 0.5 according to the previous literatures [74–77] in our experiment. For disease$d_i$$d_i$, the contribution to itself is 1 and the contribution from its children decreased as the distance between them increased. Since disease $d_i$$d_i$ lies in the most inner layer of its DAG, it is the most specific disease term whose contribution to semantic value of itself is defined as 1. Disease located on the outer layer is considered to be a more general disease, so its contribution is multiplied by the semantic contribution decay factor.

Then, through summing up all contributions from$d_i$$d_i$’s ancestors, we got the semantic value of$d_i$$d_i$. It is obvious that two diseases will get a greater similarity score when they have a larger shared part of their DAGs. Thus, we could define the $nd\times nd$$nd\times nd$ disease semantic similarity matrix 1 between disease $d_i$$d_i$ and $d_j$$d_j$ as follows:

(3) $SS1\left(d_i,d_j\right)=\frac{{\sum }_{t\in D\left(d_i\right)\cap D\left(d_j\right)}\left(D{1}_{d_i}\left(t\right)+D{1}_{d_j}\left(t\right)\right)}{DV1\left(d_i\right)+DV1\left(d_j\right)}$$SS1\left(d_i,d_j\right)=\frac{{\sum }_{t\in D\left(d_i\right)\cap D\left(d_j\right)}\left(D{1}_{d_i}\left(t\right)+D{1}_{d_j}\left(t\right)\right)}{DV1\left(d_i\right)+DV1\left(d_j\right)}$(3)

Here, the (i, j) entry of SS1 denotes the semantic similarity score between $d_i$$d_i$ and$d_j$$d_j$.

Disease semantic similarity 2

According to the existing computational model [39], we calculated the diseases semantic similarity in the second semantic similarity model.

However, it is obvious that different disease terms in the same layer of DAG($d_i$$d_i$) may also appear in some other diseases’ DAG. For example, if two diseases appear in the same layer of DAG($d_i$$d_i$) and the first disease appears in less disease DAGs than the second disease, we could conclude that the first disease is more specific than the second disease. Thus, according to the above consideration, assigning the same contribution value to these two diseases may not completely accurate since a more specific disease should have a greater contribution to the semantic value of disease$d_i$$d_i$. Therefore, the contribution of semantic value to disease $d_i$$d_i$ from the first disease should be higher than the second disease. The contribution of disease t to the semantic value of disease $d_i$$d_i$ in DAG($d_i$$d_i$) was calculated as:

(4) $D{2}_{d_i}\left(t\right)=-log\left(\frac{thenumberofDAGsincludingt}{thenumberofdiseases}\right)$$D{2}_{d_i}\left(t\right)=-log\left(\frac{thenumberofDAGsincludingt}{thenumberofdiseases}\right)$(4)

Then, we used the semantic value of disease $d_i$$d_i$ and$d_j$$d_j$, $DV2\left(d_i\right)$$DV2\left(d_i\right)$and$DV2\left(d_j\right)$$DV2\left(d_j\right)$, to compute the semantic similarity $SS2\left(d_i,d_j\right)$$SS2\left(d_i,d_j\right)$ between disease $d_i$$d_i$ and $d_j$$d_j$ based on the second model. And the semantic value DV2 was calculated by the same way as the first disease semantic similarity model by Equation (2). Further, the corresponding disease semantic similarity matrix SS2 was constructed as follows:

(5) $SS2\left(d_i,d_j\right)=\frac{{\sum }_{t\in D\left(d_i\right)\cap D\left(d_j\right)}\left(D{2}_{d_i}\left(t\right)+D{2}_{d_j}\left(t\right)\right)}{DV2\left(d_i\right)+DV2\left(d_j\right)}$$SS2\left(d_i,d_j\right)=\frac{{\sum }_{t\in D\left(d_i\right)\cap D\left(d_j\right)}\left(D{2}_{d_i}\left(t\right)+D{2}_{d_j}\left(t\right)\right)}{DV2\left(d_i\right)+DV2\left(d_j\right)}$(5)

Finally, we combined the two semantic similarity matrixes SS1 and SS2 by averaging them to obtain the disease sematic similarity matrix SS in this study.

(6) $SS\left(d_i,d_j\right)=\frac{SS1\left(d_i,d_j\right)+SS2\left(d_i,d_j\right)}{2}$$SS\left(d_i,d_j\right)=\frac{SS1\left(d_i,d_j\right)+SS2\left(d_i,d_j\right)}{2}$(6)

Gaussian interaction profile kernel similarity for miRNAs

As explained in the previous literature [78], the Gaussian interaction profile kernel similarity between miRNA $m_i$$m_i$ and miRNA $m_j$$m_j$ is constructed as follows based on the assumption that similar diseases tend to be associated with functionally similar miRNAs. The binary interaction profile vectors IP($m_i$$m_i$) and IP($m_j$$m_j$) are respectively defined as the i-th row vector and the j-th row vector of the matrix A which means whether two miRNAs are related to each disease or not. Then, we constructed the Gaussian interaction profile kernel similarity Matrix KM as follows:

(7) $KM\left(m_i,m_j\right)=e^{-{\gamma }_{m}IP\left(m_i\right)-IP{\left(m_j\right)}^2}$$KM\left(m_i,m_j\right)=e^{-{\gamma }_{m}IP\left(m_i\right)-IP{\left(m_j\right)}^2}$(7)

Here, the ${\gamma }_m$${\gamma }_m$ is used to control the kernel bandwidth which is defined by normalizing a bandwidth parameter ${{\gamma }_m}^{\prime }$${{\gamma }_m}^{\prime }$ divided by the average number of diseases associated with each miRNA. To simplify our calculation, we set the value of ${{\gamma }_m}^{\prime }$${{\gamma }_m}^{\prime }$ to 1 in our study as in the former literatures [33,78].

(8) ${\gamma }_m=\frac{{\gamma }_m^{\prime }}{\frac{1}{nm}{\sum }_{i=1}^{nm}IP{\left(m_i\right)}^2}$${\gamma }_m=\frac{{\gamma }_m^{\prime }}{\frac{1}{nm}{\sum }_{i=1}^{nm}IP{\left(m_i\right)}^2}$(8)

Gaussian interaction profile kernel similarity for diseases

Based on the assumption that functionally similar miRNAs tend to be associated with similar diseases, we could use the similar method as above to calculate the disease’s Gaussian interaction profile kernel similarity matrix KD as:

(9) $KD\left(d_i,d_j\right)=e^{-{\gamma }_{d}IP\left(d_i\right)-IP{\left(d_j\right)}^2}$$KD\left(d_i,d_j\right)=e^{-{\gamma }_{d}IP\left(d_i\right)-IP{\left(d_j\right)}^2}$(9)

where binary interaction profile vectors IP($d_i$$d_i$) and IP($d_j$$d_j$) are defined as the i-th column vector and the j-th column vector of the matrix A. The ${\gamma }_d$${\gamma }_d$ in the equation is defined by normalizing a bandwidth parameter ${{\gamma }_d}^{\prime }$${{\gamma }_d}^{\prime }$ divided by the average number of miRNAs associated with each disease. Likewise, we still set the value of ${{\gamma }_d}^{\prime }$${{\gamma }_d}^{\prime }$ to 1 in our study as in the literatures [33,78].

(10) ${\gamma }_d=\frac{{{\gamma }^{\prime }}_d}{\frac{1}{nd}{\sum }_{i=1}^{nd}IP{\left(d_i\right)}^2}$${\gamma }_d=\frac{{{\gamma }^{\prime }}_d}{\frac{1}{nd}{\sum }_{i=1}^{nd}IP{\left(d_i\right)}^2}$(10)

Integrated similarity for miRNAs and diseases

Based on the consideration that the miRNA functional similarity scores do not cover all miRNAs, we integrated the miRNA functional similarity matrix FS and the miRNA Gaussian interaction profile kernel similarity matrix KM to form a more comprehensive similarity repression. In detail, if a miRNA pair has miRNA functional similarity score, we directly considered the score as integrated similarity. On the contrary, if a miRNA pair has no miRNA functional similarity score, we utilized the Gaussian interaction profile kernel similarity score as integrated similarity. Thus, we calculated an integrated similarity matrix SM for miRNAs as:

(11) $SM\left(m_i,m_j\right)=\left\{\begin{array}{c}FS\left(m_i,m_j\right)if{m}_{i}and{m}_{j}havefunctional\\ similarityKM\left(m_i,m_j\right)otherwise\end{array}\right.$$SM\left(m_i,m_j\right)=\left\{\begin{array}{c}FS\left(m_i,m_j\right)if{m}_{i}and{m}_{j}havefunctional\\ similarityKM\left(m_i,m_j\right)otherwise\end{array}\right.$(11)

In the same way, we integrated the disease semantic similarity matrix SS and the disease Gaussian interaction profile kernel similarity matrix KD. If disease $d_i$$d_i$ and $d_j$$d_j$ have their own DAGs (i.e. these two diseases have semantic similarity), then the final integrated similarity is the average between two semantic similarity matrixes SS1 and SS2. Otherwise, the integrated disease similarity equals to the value of Gaussian interaction profile kernel similarity between them. Therefore, the integrated similarity matrix SD for disease was calculated as:

(12) $SD\left(d_i,d_j\right)=\left\{\begin{array}{c}SS\left(d_i,d_j\right)if{d}_{i}and{d}_{j}havesemanticsimilarity\\ KD\left(di,d_j\right)otherwise\end{array}\right.$$SD\left(d_i,d_j\right)=\left\{\begin{array}{c}SS\left(d_i,d_j\right)if{d}_{i}and{d}_{j}havesemanticsimilarity\\ KD\left(di,d_j\right)otherwise\end{array}\right.$(12)

Jaccard-similarity for miRNAs

The previous research [79] showed that, for recommender systems, the interaction data about customers and products could be more informative than the metadata about consumers and products. As mentioned above, the Gaussian interaction profile kernel similarity was calculated based on the known miRNA-disease association data in our study. Here we further introduced another way to represent the information of known association data.

We additionally constructed a miRNA Jaccard-similarity JM to present the similarity between two miRNAs based on their association information to diseases. The Jaccard-similarity between $m_i$$m_i$ and $m_j$$m_j$ is computed by the number of intersection divided by the number of union between sets of their associated diseases.

(13) $JM\left(m_i,m_j\right)=\frac{\left|DA\left(m_i\right)\cap DA\left(m_j\right)\right|}{\left|DA\left(m_i\right)\cup DA\left(m_j\right)\right|}$$JM\left(m_i,m_j\right)=\frac{\left|DA\left(m_i\right)\cap DA\left(m_j\right)\right|}{\left|DA\left(m_i\right)\cup DA\left(m_j\right)\right|}$(13)

where $DA\left(m_i\right)$$DA\left(m_i\right)$ denotes the set of associated diseases for miRNA$m_i$$m_i$.

Jaccard-similarity for diseases

Similar to the miRNAs, in order to present the similarity between diseases based on their association to miRNAs, we computed the Jaccard-similarity between $d_i$$d_i$ and $d_j$$d_j$ by the number of intersection divided by the number of union between the sets of their associated miRNAs. Thus, the disease Jaccard-similarity JD can be constructed as:

(14) $JD\left(d_i,d_j\right)=\frac{\left|MA\left(d_i\right)\cap MA\left(d_j\right)\right|}{\left|MA\left(d_i\right)\cup MA\left(d_j\right)\right|}$$JD\left(d_i,d_j\right)=\frac{\left|MA\left(d_i\right)\cap MA\left(d_j\right)\right|}{\left|MA\left(d_i\right)\cup MA\left(d_j\right)\right|}$(14)

where $MA\left(d_i\right)$$MA\left(d_i\right)$ denotes the set of associated miRNAs for disease$d_i$$d_i$.

Enhanced similarity-based representation for miRNAs and diseases

In order to represent both the integrated miRNA similarity and the miRNA Jaccard-similarity, we combined these two similarities together by expanding our integrated similarity matrix with the dimensionality of $nm\mathrm{x}nm$$nm\mathrm{x}nm$ to a bigger matrix of$nm\mathrm{x}2nm$$nm\mathrm{x}2nm$. The left $nm\mathrm{x}nm$$nm\mathrm{x}nm$ square matrix was the integrated similarity matrix while the right $nm\mathrm{x}nm$$nm\mathrm{x}nm$ square matrix was the Jaccard similarity matrix. Thus, we constructed a new similarity representation – miRNA enhanced similarity-based representation matrix M.

Then, similar to the miRNAs, we combined the integrated disease similarity SD and disease Jaccard-similarity JD together to obtain disease enhanced similarity-based representation matrix D. The enhanced similarity-based representation matrix M and D will be used as the feature matrices for our model.

Bipartite local models (BLMs)

BLMs [80] consider the miRNA-disease prediction problem as a link prediction problem in bipartite graphs. There are two vertex classes in the graph: one corresponds to miRNAs while the other corresponds to diseases. An edge $e_{ij}$$e_{ij}$ in the graph corresponds to the known association between miRNA $m_i$$m_i$ and disease$d_j$$d_j$. Thus, we trained two independent local models respectively based on the miRNA perspective and disease perspective to predict the likelihood score for an unknown miRNA-disease pair$e_{ij}$$e_{ij}$. Subsequently, we further aggregated these two predictions.

To be specific, we first predicted the potential association from the miRNA perspective. The prediction is based on a specific miRNA $m_i$$m_i$ and the diseases in the graph. Each disease except $d_j$$d_j$ which was being predicted was labeled as 1 or 0 depending on whether or not the $m_i$$m_i$ already had a known association with it. Then we used these labeled data to train our model to classify diseases into 1-labeled or 0-labeled classes. Subsequently, we used this model to predict the likelihood score of the investigated miRNA-disease pair$e_{ij}$$e_{ij}$. The predicted outcome for the potential association between the miRNA $m_i$$m_i$ and the disease $d_j$$d_j$ from this first model was denoted by$y_1\left(m_i,d_j\right)$$y_1\left(m_i,d_j\right)$.

Similar to the first prediction, we then predicted the potential association from the disease perspective based on a specific disease $d_j$$d_j$ and the miRNAs in the graph. However, we instead used the labeled miRNAs except $m_i$$m_i$ which are obtained according to their known association with $d_j$$d_j$ to train the model. Consequently, we used this model to predict the likelihood score of the unknown association between $m_i$$m_i$ and$d_j$$d_j$. The predicted outcome from this second model is denoted by$y_2\left(m_i,d_j\right)$$y_2\left(m_i,d_j\right)$.

Finally, we needed an aggregation function such as maximum or minimum and so on to obtain the final outcome. In this study, we chose to average the outcome from the above two models to obtain our final prediction score of the BLM for the potential association between the $m_i$$m_i$ and$d_j$$d_j$.

(15) $y\left(m_i,d_j\right)=\frac{y_1\left(m_i,d_j\right)+y_2\left(m_i,d_j\right)}{2}$$y\left(m_i,d_j\right)=\frac{y_1\left(m_i,d_j\right)+y_2\left(m_i,d_j\right)}{2}$(15)

Prediction based on weighted profiles

The BLM method has a shortcoming that it can’t handle the new diseases or miRNAs which have no known associations with miRNAs or diseases in the training data. Thus, for a new miRNA, all diseases except the predicted disease will be labeled as 0, and the local model will fail to make a prediction in this case. Similar to the new miRNA, a new disease will face the same problem that the local model can’t work with the all-0-labeled training data.

To solve this problem, when we predicted for a new miRNA$m_i$$m_i$, we used the assumption that similar miRNAs are likely to be associated with same diseases. Thus, we used the similarity between $m_i$$m_i$and other miRNAs which have known associations with $d_j$$d_j$ to calculate the likelihood score between miRNA $m_i$$m_i$ and disease$d_j$$d_j$. Therefore, the prediction score $y_1\left(m_i,d_j\right)$$y_1\left(m_i,d_j\right)$ from miRNA perspective can be computed by the weighted average of the known association between each miRNA except $m_i$$m_i$ and the disease $d_j$$d_j$ which is weighted by the similarity between $m_i$$m_i$ and other miRNAs.

(16) $y_1\left(m_i,d_j\right)=\frac{{\sum }_{m^{\prime }\in M\mathrm{\setminus }\left\{m_i\right\}}sim\left(m_i,m^{\prime }\right)i\left(m^{\prime },d_j\right)}{{\sum }_{m^{\prime }\in M\mathrm{\setminus }\left\{m_i\right\}}sim\left(m_i,m^{\prime }\right)}$$y_1\left(m_i,d_j\right)=\frac{{\sum }_{m^{\prime }\in M\mathrm{\setminus }\left\{m_i\right\}}sim\left(m_i,m^{\prime }\right)i\left(m^{\prime },d_j\right)}{{\sum }_{m^{\prime }\in M\mathrm{\setminus }\left\{m_i\right\}}sim\left(m_i,m^{\prime }\right)}$(16)

where M denotes the set of all miRNAs in the training set, $m\prime$$m\prime$means one of all other miRNAs except$m_i$$m_i$, $sim\left(m_i,m\prime \right)$$sim\left(m_i,m\prime \right)$means the similarity between miRNA $m_i$$m_i$ and $m\prime$$m\prime$ which could be obtained by the matrix M, and $i\left(m^{\prime },d_j\right)$$i\left(m^{\prime },d_j\right)$ means whether or not the miRNA $m\prime$$m\prime$ has a known association with the disease $d_j$$d_j$ which is presented in the adjacency matrix A.

$i\left(m_i,d_j\right)=\left\{\begin{array}{c}1if{m}_{i}and{d}_{j}hasaknownassociation\\ 0otherwise\end{array}\right.$$i\left(m_i,d_j\right)=\left\{\begin{array}{c}1if{m}_{i}and{d}_{j}hasaknownassociation\\ 0otherwise\end{array}\right.$

Similar to the case of the new miRNA, we predicted for a new disease $d_j$$d_j$ by the weighted average of the known association between each disease except $d_j$$d_j$ and the miRNA $m_i$$m_i$ which is weighted by the similarity between $d_j$$d_j$ and other diseases.

(17) $y_2\left(m_i,d_j\right)=\frac{{\sum }_{d^{\prime }\in D\mathrm{\setminus }\left\{d_j\right\}}sim\left(d_j,d^{\prime }\right)i\left(m_i,d^{\prime }\right)}{{\sum }_{d^{\prime }\in D\mathrm{\setminus }\left\{d_j\right\}}sim\left(d_j,d^{\prime }\right)}$$y_2\left(m_i,d_j\right)=\frac{{\sum }_{d^{\prime }\in D\mathrm{\setminus }\left\{d_j\right\}}sim\left(d_j,d^{\prime }\right)i\left(m_i,d^{\prime }\right)}{{\sum }_{d^{\prime }\in D\mathrm{\setminus }\left\{d_j\right\}}sim\left(d_j,d^{\prime }\right)}$(17)

where D denotes the set of all diseases in the training set, $d\prime$$d\prime$means one of all other diseases except$d_i$$d_i$, $sim\left(d_j,d^{\prime }\right)$$sim\left(d_j,d^{\prime }\right)$means the similarity between disease $d_j$$d_j$ and $d\prime$$d\prime$ which could be obtained by the matrix D, and $i\left(m_i,d\prime \right)$$i\left(m_i,d\prime \right)$ means whether or not the miRNA $m_i$$m_i$ has a known association with the disease $d\prime$$d\prime$ which is presented in the adjacency matrix A.

k-nearest neighbor regression with error correction (ECkNN)

BLM is a generic framework in which various regressors or classifiers can be used as local models. Bleakley et al. and Yamanishi et al. [80] used support vector machines with a domain-specific kernel in their study. In our study, we chose a hubness-aware regression model, ECkNN [81], as our local model. The algorithm of k-nearest neighbors is one of the most popular regression techniques. When using kNN to predict potential association between a miRNA and a disease from the miRNA perspective, we first determined the k-nearest neighbors of this miRNA by the ‘distance’ between them, and then use the k-nearest miRNAs to compute the predicted score of this association.

However, although the k-nearest neighbor regression has numerous advantages as described above, the algorithm has recently been found to have a drawback. The presence of ‘bad hubs’ will have negative influence on the performance of the kNN. Therefore, we used an error correction technique to alleviate the negative influence of bad hubs. Firstly, we defined the corrected label $y_c\left(x\right)$$y_c\left(x\right)$ of a training instance x as:

(18) $y_c\left(x\right)=\left\{\begin{array}{c}\frac{1}{\left|R_x\right|}\sum _{x_i\in R_x}y\left(x_i\right)if\left|R_x\right|\ge 1\\ y\left(x\right)otherwise\end{array}\right.$$y_c\left(x\right)=\left\{\begin{array}{c}\frac{1}{\left|R_x\right|}\sum _{x_i\in R_x}y\left(x_i\right)if\left|R_x\right|\ge 1\\ y\left(x\right)otherwise\end{array}\right.$(18)

where $y\left(x\right)$$y\left(x\right)$ denotes the original label of x which is uncorrected. The uncorrected original labels are directly obtained from the known miRNA-disease associations in the database. For the known miRNA-disease associations, the values of the labels were set to 1.0, and the values of the labels of other miRNA-disease associations were set to 0. $R_x$$R_x$is the set of the ‘reversed neighbors’ of x which means the set of instances whose k-nearest neighbors include x.

Then, we used the corrected labels to make prediction by the ECkNN. The predicted label $\hat{y}\left(x\prime \right)$$\hat{y}\left(x\prime \right)$ for an unlabeled instance $x\prime$$x\prime$ could be computed using the formula below:

(19) $\hat{y}\left(x^{\prime }\right)=\frac{1}{k}\sum _{x_i\in kNN\left(x^{\prime }\right)}{y}_c\left(x_i\right)$$\hat{y}\left(x^{\prime }\right)=\frac{1}{k}\sum _{x_i\in kNN\left(x^{\prime }\right)}{y}_c\left(x_i\right)$(19)

where $kNN\left(x^{\prime }\right)$$kNN\left(x^{\prime }\right)$means the set of the k-nearest neighbors of$x\prime$$x\prime$. And the k was set to 100 in our experiment according to the literature [81] that the mean absolute error gradually decreased and eventually converged as the k increased in their experiment.

Blharmda

In this study, we developed BLHARMDA to predict potential miRNA-disease associations. Our model could be divided into three steps. The whole process was illustrated in Figure 3. The first step of this model is preparing the data for the following prediction. Then, the second step of our model is using the data we obtained in the first step to predict the potential associations by BLM or weighted profile method. We only used weighted profile method on the new miRNAs or diseases without any known related diseases or miRNAs. For other miRNAs or diseases with known related diseases or miRNAs, we made predictions based on the BLM method with a hubness-aware regression, ECkNN regression, as its local model.

Figure 3. Flowchart of potential miRNA-disease association prediction based on the computational model of BLHARMDA: 1) data preparation, where enhanced similarity representation for miRNAs and diseases were constructed in this step; 2) Training the BLM with ECkNN as the local model and making predictions for the miRNAs or diseases with known associations by the BLM we trained. Then, using the weighted profile method to make prediction for new miRNAs and diseases.

In the BLM method, we predicted the associations by two models based on the miRNA perspective and disease perspective respectively. Subsequently, through using the aggregating function of arithmetical average on the two prediction outcomes above, we obtained our final prediction score ${\hat{y}}_k\left(m_i,d_j\right)$${\hat{y}}_k\left(m_i,d_j\right)$ of an association by BLM with ECkNN as the local model.

In the weighted profile method for new miRNAs and diseases, we also predicted the associations by two models from the miRNA perspective and disease perspective respectively. Still, we only used the selected features of matrix M and D to obtain our two prediction outcomes. Then, our final prediction score ${\hat{y}}_w\left(m_i,d_j\right)$${\hat{y}}_w\left(m_i,d_j\right)$ of this association can be obtained by averaging these two outcomes.

Thereafter, by integrating the prediction scores from the BLM and the prediction scores of new miRNAs and diseases obtained by the weighted profile method, we could get a global prediction of the associations between all miRNAs and diseases. Thus, the prediction outcome $\hat{y}\left(m_i,d_j\right)$$\hat{y}\left(m_i,d_j\right)$ can be computed as follows:

(20) $\hat{\mathrm{y}}\left({\mathrm{m}}_{\mathrm{i}},{\mathrm{d}}_{\mathrm{j}}\right)=\left\{\begin{array}{c}{\hat{\mathrm{y}}}_{\mathrm{w}}\left({\mathrm{m}}_{\mathrm{i}},{\mathrm{d}}_{\mathrm{j}}\right)if{m}_{i}or{d}_{j}isnewmiRNAordisease\\ {\hat{\mathrm{y}}}_{\mathrm{k}}\left({\mathrm{m}}_{\mathrm{i}},{\mathrm{d}}_{\mathrm{j}}\right)if{m}_{i}or{d}_{j}is{n}^{\prime }tnewmiRNAordisease\end{array}\right.$$\hat{\mathrm{y}}\left({\mathrm{m}}_{\mathrm{i}},{\mathrm{d}}_{\mathrm{j}}\right)=\left\{\begin{array}{c}{\hat{\mathrm{y}}}_{\mathrm{w}}\left({\mathrm{m}}_{\mathrm{i}},{\mathrm{d}}_{\mathrm{j}}\right)if{m}_{i}or{d}_{j}isnewmiRNAordisease\\ {\hat{\mathrm{y}}}_{\mathrm{k}}\left({\mathrm{m}}_{\mathrm{i}},{\mathrm{d}}_{\mathrm{j}}\right)if{m}_{i}or{d}_{j}is{n}^{\prime }tnewmiRNAordisease\end{array}\right.$(20)

In order to promote further research of biologists, we published our source code in https://github.com/BLHARMDA/Bipartite-Local-Models-and-hubness-aware-regression.

Disclosure statement

No potential conflict of interest was reported by the authors.

Supplementary material

Supplemental data for this article can be accessed here.

Additional information

Funding

XC was supported by National Natural Science Foundation of China under Grant No. 61772531.