MGOGP: a gene module-based heuristic algorithm for cancer-related gene prioritization

Discovering cancer-related genes has profound applications in modelling, diagnosis, therapeutic intervention, and in helping researchers get clues on which genes to explore [1–3]. Computational approaches are preferred due to their high efficiency and low cost [4, 5]. Many computational methods have been proposed, including: a) gene-based function similarity measure methods [6–9]; b) biological interaction network-based methods [10–14], and c) methods based on multiple datasets fusion [15–17]. Methods of the first kind based on the hypothesis that phenotypically similar diseases are caused by functionally related genes. Based on this hypothesis, many methods prioritize genes by computing similarity scores between the candidate genes and the known disease genes. For example, ToppGene [6] ranks genes based on similarity scores of each annotation of each candidate genes by comparing enriched terms in a given set of training genes. Endeavour [8] prioritizes candidate genes by similarity values between candidate genes and seed genes, by integrating more than six types of genomic datasets from over a dozen data sources. Methods of the second kind prioritize genes using the guilt-by-association principle, which means genes interacting with known disease genes are more likely disease-related genes. For instance, PINTA [10] prioritizes candidate genes by utilizing an underlying global protein interaction network. Other methods rank candidate genes by exploiting either local or global network information [2]. Methods of the last kind incorporate datasets such as gene expression, biomedical literature, gene ontology, and PPIs together for gene prioritization. For example, ProphNet [17] integrates information of different types of biological entities in a number of heterogeneous data networks. Taking all these methods into consideration, they are either gene-centric or network-centric. However, gene module as a basic functional unit of genes has seldom been considered.

Gene module can be defined as a protein complex, a pathway, a sub-network of protein interactions. Module detection has long been studied and many useful algorithms have been proposed, such as [18–21]. Although different methods have different module detection strategies, most of them rely on PPIs network. PPIs network suffers from drawbacks as highlighted in [22]. Firstly, the PPI network is incomplete, which only covers the interactions of well-researched proteins. For instance, of the 20,502 genes in the gene expression matrix downloaded from The Cancer Genome Atlas (TCGA), only 9078 (44.2%) and 2761 (13.4%) genes are included in Human Protein Reference Database (HPRD) [23] and Database of Interacting Proteins (DIP) [24] PPIs networks respectively. As a result, detected modules are incomplete and their accuracy are limited. Secondly, protein interactions in PPIs network suffer from high false positive and negative rates, modules discovered from such PPI data also suffer from high false rates. All these inherent limitations affect the coverage and accuracy of the inferred modules.

Nowadays, numerous public databases of protein and gene annotation information are available, such as Entrez Gene [25], Ensembl [26], PIR iProClass [27], GeneCards [28], KEGG [29], Gene Ontology Consortium [30], DAVID [31], GSEA [32] and UniProt [33]. For instance, DAVID [31] contains information on over 1.5 million genes from more than 65,000 species, with annotation types, including sequence features, protein domain information, pathway maps, enzyme substrates and reaction, protein-protein interaction data and disease associations. Gene Ontology Consortium describes the functions of specific genes, using terms known as GO (Gene Ontology). KEGG map genes to pathways while GSEA provides functional gene groups collected from BioCarta genes sets, KEGG gene sets and Reactome gene sets. With these annotation information, we can easily group genes into functional modules.

Complex diseases, especially cancer are caused by the dysfunction of groups of genes and/or gene interactions rather than the mutations of individual genes. Detecting and prioritizing cancer-related genes from the perspective of gene module is promising. Although some useful work has been conducted [34, 35], the results are still far from being satisfactory. In this study, we take the importance of not only genes but also their affiliated modules into consideration, and prioritizing genes in a heuristic way. We measure module importance by the number of differential genes within the module and the number of differential correlations between the module genes. Besides, the number of known cancer-related genes in the module is also considered. We measure the gene importance by three aspects information: a), gene’s differential expression value, b), the number of differential correlations between the gene and all other module gene. c), the fuzzy measure based similarity values between the gene and all known cancer-related genes (if exist) within the module. The global rank of all genes is obtained by utilizing a rank fusion strategy.

As shown in Fig.1, MGOGP takes gene expression datasets, gene modules, known disease genes and gene ontology annotation information [36] as input, and the ranked genes as output. The main parts including: module importance measure, module-specific gene importance measure, module rank and module-specific gene prioritization, and global cancer-related gene prioritization. Figure 2 schematically illustrates these steps in detail.

First, obtain functional gene modules; then get the global ranking of all modules and the local ranking of all module-specific genes based on their importance; finally, the rank fusion algorithm further gives all genes a global rank.

Global gene ranking

Most genes deploy their functions in the context of sophisticated functional modules [45, 46]. Therefore, the global rank of a gene need be decided by its own importance and the importance of its affiliated module. As in [34], a rank fusion strategy is used to fuse the local rank of genes in each module into a global rank.

The rank fusion strategy is a recursive process. It decides the rank of the nth gene based on all the top-ranked n − 1 genes. We define i as the number of genes having already obtained their global ranking in the recursive process of rank fusion, m(i, j) as the number of top i genes located in the module j after having determined the top i genes. t(i, j) as the expectation of the number of top i genes located in the module j. e(i, j) as the expectation of probability that the i + 1 globally ranked genes come from the module j. We use the module importance value p(m _j ) as the probability of a disease-related gene comes from it. The relationship between i, m(i, j), t(i, j) and p(m _j ) is shown in Eq. 20:

$$ {\displaystyle \begin{array}{l}t\left(i,j\right)= ip\left({m}_j\right)\\ {}e\left(i,j\right)=t\left(i+1,j\right)-m\left(i,j\right)\end{array}} $$

(20)

Initially, the first ranked gene in the module with highest importance value is chosen as the top 1 gene in the gene’s global rank, because all genes in each module have been ranked from big to small according to their importance value. Let i as the number of genes having obtained their global ranking, to decide the i + 1 ranked gene, we need to find the module with the biggest e(i, j) value, because e(i, j) indicates the expectation of probability that the i + 1 globally ranked genes from module j. So the genes ranked m(i, j) + 1 in the module j will be chosen as the top i + 1 ranked gene, because in the module j, top m(i, j) genes has obtained the global ranking. Repeat the process until all genes get ranked. As shown in Fig. 3 (in Additional file 1).

Both raw count and normalized gene expression datasets are downloaded from TCGA (http://cancergenome.nih.gov/) [47], which include expression values of 20,503 genes across 102 normal samples and 779 tumor samples. Besides, gene expression datasets of Prostate adenocarcinoma containing 483 tumor samples and 51 normal samples are also downloaded from TCGA. Four thousand seven hundred twenty-six gene modules are downloaded from the website of GSEA (in Additional file 2).

Next, we use MGOGP for genome-wide breast cancer gene prioritization. We use 328 breast disease-related genes downloaded from SNP4Disease (http://snp4disease.mpibn.mpg.de/result.php) as seed genes (see Additional file 3). Ten well-known breast cancer-related genes (shown in Table 4, which are not contained in the 328 genes) are used to validate the effectiveness of our method. All GSEA gene modules are pre-processed by removing all the genes which do not have gene expression information (the final module list is supplied in Additional file 4). The result is shown in Fig. 4.

Gene ID	Gene symbol	Gene name
672	BRCA1	Breast Cancer 1, Early Onset
675	BRCA2	Breast Cancer 2, Early Onset
7157	TP53	Tumor Protein P53
5728	PTEN	Phosphatase And Tensin Homolog
841	CASP8	Caspase 8, Apoptosis-Related Cysteine Peptidase
2263	FGFR2	Fibroblast Growth Factor Receptor 2
4214	MAP3K1	Mitogen-Activated Protein Kinase Kinase Kinase 1, E3 Ubiquitin Protein Ligas
11200	CHEK2	Checkpoint Kinase 2
472	ATM	ATM Serine/Threonine Kinase
83990	BRIP1	BRCA1 Interacting Protein C-Terminal Helicase 1

As can be seen from Table 5, many top-ranked modules are included in well-known breast cancer pathways, such as PI3K/AKT [48] pathway and VEGF ligand-receptor pathway. The VEGF family of ligands and receptors are intimately involved in tumor angiogenesis, lymphangiogenesis, and metastasis [49]. More importantly, of the 100 genes in the top 10 ranked modules, 20 of them are contained in the KEGG breast cancer pathway (hsa05224), which is an indication of the good performance of MGOGP for cancer gene prioritization.

Next, we validate the performance of MGOGP by comparing the gene prioritization results with results obtained by methods: Endeavour [8], GeneFriends [50], PINTA [10], TOPPGene [6] and TOPNet [13]. All the methods use the same datasets and under their default parameter settings. The results are shown in Fig. 5. Brief descriptions of these methods are provided in Additional file 6. Core sourcecode of MGOGP is provided in Additional file 7. Other source codes are available from the corresponding author on reasonable request.

In Fig. 5, we count the number of breast cancer-related genes in the gene prioritization results. As is shown in Fig. 5, MGOGP outperforms other methods in detecting cancer-related genes. We use all the 328 breast disease related genes as known disease gene (Endeavour and GeneFriends used the same gene sets) and count the number of known disease genes appear in top 100–1000 prioritization results.

To do comparison more rigorously, we further compare MGOGP to Endeavour, TOPNet and TOPPGene. Each time we randomly select 100, 150 and 200 different known disease genes from the 328 breast disease-related genes for known disease genes and others are left for test (each kind of selection repeat 100 times). We count the average number of test genes appear in Top 200 gene prioritization results. Results are shown in Fig. 6.

In Fig. 7, the number of Known Disease Gene is the number of genes supplied for training each method that fall within the top 10. For example, in Table 6, PTEN, VEGFB, and MCM2 are three genes fall within the top 10 of the gene ranking result, so the number of Known Disease Gene of MGOGP in Fig. 7 is 3. For each gene within the top 10 gene ranking results of each method, we search the number of articles in PubMed mention the association between the gene and breast cancer. We count the number of genes has more than 10 PubMed article reference. As shown in Fig. 7, genes detected by MGOGP have more article supports than other methods.

Results of omics experiments commonly consist of a large set of genes, while researchers usually only care about the behaviour of several genes. In this paper, a heuristic algorithm is proposed for prioritizing disease-associated genes by utilizing gene ontology annotation information and known disease-related genes as heuristic information. Different from existing methods, we propose to rank genes considering the importance of both individual genes and their affiliated modules, and utilize Gene Ontology (GO) based fuzzy measure value as well as known disease genes as heuristics, and use rank fusion strategy to obtain the global gene prioritization. Results show that MGOGP outperforms many other methods in cancer-related gene prioritization.

Different from other module-based gene prioritization methods, where modules are detected by partitioning the network using the network clustering methods, we obtain modules through gene function annotation, that is, functionally related genes are grouped into the same modules. Because gene interaction networks often suffer from the problems of high rates of false positive/negative interactions, and modules detected by network clustering algorithms often have limited accuracy, so our method is more advanced. One important difference between modules used in this study and modules detected through network partition is that no edges in our module. Instead, we use statistical methods detecting differential correlations between genes within a module, which could help avoid the preference of genes or modules that are well-researched (because currently obtained network is far from complete, the number of interactions among well-researched genes may be much more than that of newly discovered genes).

Different from module-based methods in [34], MGOGP ranks modules considering three aspects of information: module-specific gene importance, differential correlations, and importance of the module itself. In [34], the author considers the importance of a module by considering only the number of disease genes and the size of the module, which may bias toward big modules. Furthermore, gene as the major component of the module whose importance is not considered when measuring the importance of a module in [34]. While in our method, when measuring the importance of a module, we consider: the importance of the module itself, the importance of module contained genes as well as differential correlations within the module, which are the main improvements of our method.

Compared with other non-module-based prioritization methods, our algorithm also has obvious advantages. First, it is easier to find the potential pathogenic genes that cause the disease from the point of view of gene modules. Second, it takes cross-validation strategy which could guarantee the stability of the recognition results. And our method works with heuristic information which could effectively avoid the blindness of the search.

By applying MGOGP on different datasets, we demonstrate that MGOGP performs better than previous gene or network-centric methods in terms of potential disease-related genes prediction. Firstly, the performance of MGOGP is validated by comparing it with three module based cancer-related gene prioritization methods. Results show that all test genes are ranked on average within top10% of all the candidate genes. According to our results, many top-ranked modules are included in well-known cancer pathways, and top-ranked genes have more supporting PubMed articles. All of the results show that our methods perform better than the state of the art methods.

Prioritization methods are useful for assisting scientists at early research stages, and to formulate novel hypotheses of interest. In the future, one of our main goals is to see how our method behaves in other prioritization problems when using different entities and sources of data sets not covered in this study. Furthermore, we plan to study in more detail the quality of the datasets and their influence on algorithm performance, and design new methods to try to improve the results. As we all know that the methods become more mature the results will become increasingly accurate and more biologically meaningful.