A Knowledge-Based Machine Learning Approach to Gene Prioritisation in Amyotrophic Lateral Sclerosis

Amyotrophic lateral sclerosis is a neurodegenerative disease of the upper and lower motor neurons resulting in death from neuromuscular respiratory failure, typically within two to five years of first symptoms. Several rare disruptive gene variants have been associated with ALS and are responsible for about 15% of all cases. Although our knowledge of the genetic landscape of this disease is improving, it remains limited. Machine learning models trained on the available protein–protein interaction and phenotype-genotype association data can use our current knowledge of the disease genetics for the prediction of novel candidate genes. Here, we describe a knowledge-based machine learning method for this purpose. We trained our model on protein–protein interaction data from IntAct, gene function annotation from Gene Ontology, and known disease-gene associations from DisGeNet. Using several sets of known ALS genes from public databases and a manual review as input, we generated a list of new candidate genes for each input set. We investigated the relevance of the predicted genes in ALS by using the available summary statistics from the largest ALS genome-wide association study and by performing functional and phenotype enrichment analysis. The predicted sets were enriched for genes associated with other neurodegenerative diseases known to overlap with ALS genetically and phenotypically, as well as for biological processes associated with the disease. Moreover, using ALS genes from ClinVar and our manual review as input, the predicted sets were enriched for ALS-associated genes (ClinVar p = 0.038 and manual review p = 0.060) when used for gene prioritisation in a genome-wide association study.


Introduction
Amyotrophic lateral sclerosis (ALS) is a neurodegenerative disease of the upper and lower motor neurons resulting in death from neuromuscular respiratory failure, typically within 2-5 years of the first symptoms [1]. Genetic factors are an important cause of ALS, and over 130 genes [2] have been linked to the disease, although only a subset have strong and replicated supporting evidence [3]. Rare disruptive variants are responsible for about two-thirds of familial cases [4] and 10-15% of the remaining 90-95% of patients who do not report a family history of ALS [4][5][6]. Therefore, for most

Definition of ALS-Linked Genes
In total, five sets of ALS-linked genes were used in this study: (1) DisGeNet-101 genes [28], (2) tALSoD-126 genes [2,29], (3) ClinVar-44 genes [30], (4) a manually-curated list-40 genes [3] and (5) the union of all other sets-199 genes. DisGeNet contains several specific subtypes of ALS as UMLS identifiers. For this study, we chose to merge these subtypes to prevent the presence of very similar diseases making the learning unrealistically easy. Gene lists are available as Supplementary File 1. Note that while multiple definitions of ALS-linked genes were used, all other disease-gene associations were from DisGeNet.

Machine Learning
We used our recently published method for knowledge graph completion (available at https: //github.com/KHP-Informatics/ADR-graph) (version 0.1) [24]. In brief, the method takes as input a graph of known data related to the prediction task, in this case, gene-disease links, gene functions, and others, and returns a list of predicted edges missing from that graph. These relationships are predicted to exist based on similar patterns found in the graph. For example, to predict new genes linked to ALS, the algorithm compares all genes known to be linked to ALS to all other genes and builds a predictive profile based on a weighted combination of existing relationships in the graph.
Every gene is then scored for its similarity to this profile. Predictions are made by applying a threshold to this similarity score, with all genes above the threshold predicted as candidate ALS-linked genes. The optimum weighting and score threshold are learned from the known set of linked genes [24].

Cross-Validation
The goal of the machine learning algorithm is to predict new genes linked to ALS based on facts about the genes known to be associated. In the cross-validation, this task is simulated by deleting a proportion of the known genes and training the algorithm on the remaining data before testing whether it can predict the link to the missing genes correctly. Unlike a typical machine learning cross-validation, where the training and test sets are entirely separate, here, the algorithm was trained on the full set of genes, but the link to ALS of the test set was deleted from the graph. This means that during training, the "missing" linked genes were labelled as not linked to ALS. Therefore, the job of the algorithm was to replace this missing link in the graph, making the task harder as genes that were linked to ALS were included during training as true negatives. This is an exact simulation of the use-case.

Gene Set Function and Phenotype Enrichment Analyses
Gene ontology term overrepresentation was analysed using Pantherdb [31] version 15.0. All tests were performed using the "Panther GO slim for biological process" and the Fisher s exact test for overrepresentation with the Bonferroni correction for multiple testing. Enrichment vs. the background list with p-value < 0.05 after the Bonferroni correction was considered significant. Similar biological processes were grouped with Revigo [32]. For the predicted ALS-linked gene sets, the genes in the knowledge graph that were not linked to ALS according to the relevant definition were used as the background gene list (i.e., the set of all genes that could have been predicted).
We used the Enrichr webserver [33] to test our predicted gene lists for the enrichment of genes linked to human diseases in the OMIM database [34,35]. An adjusted p-value < 0.05 was considered significant.

GWAS Validation
To validate the association between the candidate genes proposed by our method and ALS, we used the largest ALS GWAS for which summary statistics were publicly available at the time [36]. It included 80,610 individuals (20,806 cases and 59,804 controls) of European ancestry. For each set of genes, only SNPs mapped onto those genes were used for validation. The single-test p-value threshold to indicate statistical significance was 0.05. When multiple tests were performed, we used the Bonferroni correction [37] based on the number of genes being tested. All annotations, genomic positions, and variants refer to the reference human genome hg19/GRCh37. The 1000 Genomes Project Phase 3 reference panel was used to compute the r 2 and MAF [38]. The gene-based and gene-set association analyses were performed with Magma v1.07b [39]. In the gene-based and gene-set analyses, the associations of all SNPs within each gene or gene-set were tested simultaneously.

Definition of Known ALS-Linked Genes
Given that there is no consensus on what genes increase the risk of ALS, we chose to separately model lists of genes from a range of sources (Supplementary File 1). These included a specialized public database of ALS gene mutations in humans, the ALS Online Database (ALSoD), a general database of clinical relevant mutations, ClinVar, a general database of all disease-associated (in the broadest sense) genes, DisGeNet, our manual review of ALS-linked genes curated by scientists working in the field, and the union of these four lists. The overlap between these lists is shown in Figure 1. Overlap between lists of ALS-linked genes from different sources. Note that areas are not to scale. There are no genes that are listed in both ClinVar and DisGeNet but not found in either ALSoD or the manual review. In total, there are 199 unique genes linked to ALS across all lists.

Model Training and Cross-Validation
We trained a model for each of the sets of ALS-linked genes separately and used 5-fold crossvalidation to estimate the performance of each predictive model ( Table 1). The purpose of the crossvalidation was to estimate how well a model trained on all known ALS-linked genes would perform in predicting novel associations, by testing whether it would have predicted some of the known genes based on the rest. Table 1. Cross-validation performance for each definition of ALS-linked genes. Precision, recall, and fold-change enrichment are given as the mean (standard deviation) for five folds. The numbers of genes are minimum-maximum. The fold-change enrichment is based on the expected precision of random guessing. Each of the five folds is considered significantly enriched if the p-value for predicting at least as many genes correctly is < 0.05 under the hypergeometric distribution.

Model Training and Cross-Validation
We trained a model for each of the sets of ALS-linked genes separately and used 5-fold cross-validation to estimate the performance of each predictive model ( Table 1). The purpose of the cross-validation was to estimate how well a model trained on all known ALS-linked genes would perform in predicting novel associations, by testing whether it would have predicted some of the known genes based on the rest. Table 1. Cross-validation performance for each definition of ALS-linked genes. Precision, recall, and fold-change enrichment are given as the mean (standard deviation) for five folds. The numbers of genes are minimum-maximum. The fold-change enrichment is based on the expected precision of random guessing. Each of the five folds is considered significantly enriched if the p-value for predicting at least as many genes correctly is < 0.05 under the hypergeometric distribution. For all lists, the training precision was low, which was expected for this task as the validation genes remained in the graph and were considered true negatives during training. Ideally, the model should have labelled these as positives (which would be considered "false positives" by standard evaluation metrics). These "false positives" were taken as a prediction that we subsequently validated. The performance of the DisGeNet list was notable. Although the precision was relatively high, the standard deviation was also much greater. This indicated that performance was very sensitive to randomisation of the folds, possibly due to a less uniform set of genes in terms of the available data. The model also achieved much lower recall for DisGeNet (only 23% in training and 9% in validation), and only 2/5 training folds were significantly enriched for validation. The lower recall of the DisGeNet list could be due to DisGeNet being a broad database, with no manual curation that also includes genes linked to ALS only in model organisms. Models based on all other lists significantly outperformed random guessing in all folds, indicating the model has learned discriminative features of the provided ALS-linked genes.

Prediction of New ALS-Linked Genes
For each set of ALS-linked genes, we trained a model on the full list and generated new predictions based on the learned profile. Taking the full set of new ALS-linked genes predicted from each list gave 45 (Manual list), 176 (DisGeNet), 192 (ClinVar), 327 (ALSoD), and 575 (union) predicted linked genes ( Figure 2). The genes predicted by the model trained on the manually curated list are shown in Table 2. The full set of predicted genes from all models are available in Supplementary File 2.
As shown in Figure 2, there were no genes predicted by all definitions of known ALS-linked genes. The main reason for this is that the input lists are not identical and frequently predict members of other lists. For example, a prediction made from ALSoD cannot be predicted from DisGeNet if it is already linked to ALS in the DisGeNet list. Moreover, 91 genes were uniquely predicted by the union set, meaning that this list did not simply produce the union of the outputs of all other lists.
In contrast, 15 genes were predicted from at least one other list but not the union list and which were not part of any list themselves. Although there are 73 genes outside of the union predictions in Figure 2, many were linked to ALS in other lists. In other words, combining all lists, and then training the model is not the same as training a model on each list and then combining the results. This is likely due to the impact of combining lists of genes on the enrichment test in the feature selection step of model training. Figure 3 shows the number of predictions from each list that were already considered to be linked to ALS in at least one other list. As shown in Figure 2, there were no genes predicted by all definitions of known ALS-linked genes. The main reason for this is that the input lists are not identical and frequently predict members of other lists. For example, a prediction made from ALSoD cannot be predicted from DisGeNet if it is already linked to ALS in the DisGeNet list. Moreover, 91 genes were uniquely predicted by the union set, meaning that this list did not simply produce the union of the outputs of all other lists.
In contrast, 15 genes were predicted from at least one other list but not the union list and which were not part of any list themselves. Although there are 73 genes outside of the union predictions in

Functional and Phenotypic Enrichment Analyses of Predicted ALS-Linked Genes
Given that our trained models were significantly enriched for ALS-linked genes in a simulation of this task (the cross-validation), the new predictions from the full models, using the whole lists as training datasets, might contain new ALS-linked genes. The predicted sets of genes overlapped to varying degrees and could be functionally similar. We used GO term overrepresentation applied to each set of predicted genes to investigate which biological processes they were enriched for. A complete list of significantly enriched processes and corresponding GO terms is available in Supplementary File 3. All predicted sets were enriched for genes involved in biological processes associated with ALS. For example, considering the 20 processes with the highest enrichment fold for practicality (Table 3), several processes related to angiogenesis [40], lipid metabolism [41], mitochondria activity [42], protein kinase activity [43], superoxide metabolism [44,45], vesicletrafficking [46], neurotransmitter regulation [47], and behaviour [48] were widely shared across the sets.  The 45 genes predicted from the Manual list of ALS-linked genes were all contained in at least one other list. All 45 were in DisGeNet, 12 were in ALSoD, and seven were in ClinVar.

Functional and Phenotypic Enrichment Analyses of Predicted ALS-Linked Genes
Given that our trained models were significantly enriched for ALS-linked genes in a simulation of this task (the cross-validation), the new predictions from the full models, using the whole lists as training datasets, might contain new ALS-linked genes. The predicted sets of genes overlapped to varying degrees and could be functionally similar. We used GO term overrepresentation applied to each set of predicted genes to investigate which biological processes they were enriched for. A complete list of significantly enriched processes and corresponding GO terms is available in Supplementary File 3. All predicted sets were enriched for genes involved in biological processes associated with ALS. For example, considering the 20 processes with the highest enrichment fold for practicality (Table 3), several processes related to angiogenesis [40], lipid metabolism [41], mitochondria activity [42], protein kinase activity [43], superoxide metabolism [44,45], vesicle-trafficking [46], neurotransmitter regulation [47], and behaviour [48] were widely shared across the sets.
We also investigated whether the sets of predicted ALS-linked genes were enriched for genes associated with other related neurodegenerative diseases. In this test, we used Enrichr with the disease-gene sets from the OMIM database (Table 4). Several neuromuscular diseases, e.g., Dystonia and Cardiomyopathy, and neurodegenerative disorders, e.g., Charcot-Marie-Tooth disease (CMTD), Parkinson s disease (PD), Frontotemporal dementia (FTD), Schizophrenia (SCZ) and Alzheimer's Disease (AD), were among the top 5 significant hits (Table 4). For some of these diseases, e.g., FTD, SCZ, PD, and CMTD, the genetic and phenotypic overlap with ALS has been shown previously [49][50][51][52][53][54][55]. Interestingly, type 2 diabetes was among the top 5 DisGeNet hits and has been proposed as an ALS risk factor and phenotypic modifier [56,57]. The predictions from the Manual list were enriched for genes already linked to ALS in the OMIM database, again highlighting the importance of using different gene lists. Table 3. Top 20 most enriched Revigo-grouped biological processes (ranked by enrichment fold) for each set of predicted ALS-linked genes. Enrichment was calculated using Pantherdb. All terms are significant at p-value < 0.05 level following the Bonferroni correction.

Validation of Predicted ALS-Linked Genes in GWAS Data
In GWAS, multiple testing correction approaches, such as the Bonferroni method, are used to protect from Type I errors. However, at a genome-wide scale, the price to pay for a minimal Type I error is the limitation of the statistical power provided by the sample. One way to overcome this issue is to restrict the number of tests to a limited number of candidate genes, such that a less stringent p-value threshold can be used. To test the relevance of the predicted genes in ALS further, and explore the feasibility of using our method to prioritise GWAS hits, we used our predicted gene sets to identify new ALS genes based on a previously published ALS GWAS. This study involved over 80,000 individuals of European ancestry, about a quarter of whom were people with ALS and the remainder unaffected controls.
Using the publicly available summary statistics, we performed a Magma gene-based analysis that combined all SNPs within each gene to test their association with ALS simultaneously. The~10,000,000 SNPs from the GWAS were mapped onto 18,067 protein-coding genes. For each list, a p-value smaller than 0.05/N was used to indicate a significant association. N was the number of genes in the set that were effectively tested in the Magma gene analysis. At least one new prediction was validated for each model by this method (Table 5). However, when we tested the likelihood of our findings occurring by chance, only the ClinVar model (p = 0.038) passed the significance threshold, while the Manual model was close (p = 0.060). Table 5. Validation of predicted ALS-linked genes. *ZFP91-CNTF is a read-through transcript of both ZPF91 and CNTF. CNTF was predicted to be ALS-linked by the model. The p-value and number of validated genes shown in brackets indicate results if the ZFP91-CNTF transcript is not considered. ATXN3 was predicted and validated for three models and WNT7A for two models. For the ZFP91-CNTF locus, the whole co-transcript, including both ZFP91 and CNTF was used instead of the two individual genes, as the SNPs within these two adjacent genes were in strong linkage disequilibrium (Supplementary Figure S1). Note that even if more than one model predicted the same gene, it would not have necessarily been validated by all of them as the p-value threshold depended on the total number of genes predicted by each model. The Magma gene-set analysis allowed us to test all the genes within each list simultaneously. Similar to the single-gene analysis, the Manual and ClinVar models were the ones with the lowest p-values (p = 0.057 and p = 0.065, respectively).

Discussion
In this study, we trained machine learning models to predict new ALS genes based on multiple available definitions of the currently known ALS genes. Overall, all models performed well in simulations and were able to predict genes involved in ALS-related diseases and biological processes. Only one of the five models (based on ClinVar) produced significant results (p = 0.038) when validated in human GWAS data.
There are several databases of disease-gene associations, as well as several literature reviews, but there is no consensus list of all currently known ALS-linked genes. Partly, this is because ALS is a complex, rare, and rapidly progressive disease whose study involves elaborate designs that are difficult to implement and replicate; partly because such lists and databases utilise different definitions of ALS-linked genes according to their aim. For example, we used a manually curated list of genes containing only genes with strong and replicated evidence of association with ALS. In contrast, large non-specialist databases might consider weaker evidence such as non-replicated differential expression in a transgenic animal model as sufficient evidence. Neither of these examples is necessarily correct or incorrect, but the choice of genes is very likely to impact the results of a predictive study such as this, where predictions are based on similarity to the genes that we use to define ALS. Our choice to build four models based on four different definitions of ALS genes and a fifth model based on the union of the four gene lists reflects such considerations.
We trained our models using only known disease-associated genes, protein-protein interactions, and gene functional annotations. We anticipate that the inclusion of additional data could significantly improve the performance of the models. However, care must be taken to ensure that additional data types do not inadvertently leak information. By this, we mean that the indirect relationships between nodes in the graph can imply information we intended to remove. For example, we considered the inclusion of drug-disease and drug-target data in the model, but as ALS is the only disease linked to Riluzole, this would link all ALS genes indirectly to Riluzole [58], making the edges deleted in the cross-validation, trivial to predict. Furthermore, some types of information are not independent. For example, a drug treatment might only be developed based on a known ALS gene, making the cross-validation challenging to evaluate fairly. In effect, this means the time series of events must be known so that it can be used to filter out subsequently derived information and exclude it from the graph, but this information is not trivial to access at scale for many data types.
Beyond the specific genes, our results also have important implications for the use of genetic databases in machine learning in general. We have shown that when used to train the models, the different definitions of known ALS genes derived from different data sources produced very different results, such that simply combining these definitions produced further distinct predictions.
The key implication of this result is that the outputs of separate studies can only be compared with considerable caution, if at all, and future work should be very transparent as to the definitions and sources of data used.
The major challenge for the use of this knowledge-based machine learning approach to predict new ALS genes is the validation of the predicted genes. We attempted several approaches to test the performance of the models and assess the relevance of the predictions in ALS. First, the cross-validation test showed that by using a subset of the genes in the lists, the models were able to predict a significant part of the missing data in most cases. The model cross-validation performance was higher for the Manual list, matching the expectation that the genes from a manually curated and conservative source might represent a more homogeneous and reliable description of the disease genetics. Secondly, we assessed which biological processes were over-represented in our predictions. All the predicted sets of genes were enriched for biological processes known to be affected by ALS including processes related to angiogenesis [40], lipid metabolism [41], mitochondria activity [42], protein kinase activity [43], superoxide metabolism [44,45], vesicle-trafficking [46] and neurotransmitter regulation [47]. Interestingly, behaviour related processes [48] were over-represented in the DisGeNet, ALSoD, and Union predictions. Behavioural changes are frequent in ALS patients [48] and closely related to the cognitive decline observed in a significant proportion of patients [59,60]. Thirdly, we performed a phenotype enrichment analysis to test whether our predictions were significantly enriched for genes linked to other diseases. This test showed that all models predicted a significant number of genes that are involved with other neurodegenerative and neuromuscular diseases whose genetic overlap with ALS is either known or for which there is increasing recognition. These included FTD, SCZ, PD, and CMTD. Again, the Manual predictions were the most conservative as the only other enriched disease was FTD whose genetic overlap with ALS is established, and the two phenotypes often coexist in patients [61].
Finally, we used our predictions in a candidate gene approach to discover new genes associated with ALS using a Magma gene analysis based on the publicly available summary statistics of the most recent ALS European GWAS. We were able to validate several associations, including some known candidate genes such as SCFD1 [62] and UNC13A [62]. However, when we considered the number of validated genes per model, only for ClinVar was this significantly larger than what was expected by chance (3 genes validated out of 170, p = 0.038). For the Manual set, this was close to significance (2 genes out of 41, p = 0.060). The Manual and ClinVar sets were also the only ones close to the significance threshold in the gene-set analysis that evaluated the association of all genes in each set and the disease simultaneously (p = 0.057 and p = 0.065 respectively). It is important to acknowledge that although the SNP based GWAS offered a possibility to validate our predictions, and machine learning methods such as ours are often used for variant prioritisation in GWAS, SNPs represent a limited amount of ALS heritability (~8%) [62]. Therefore, we expected a limited number of genes to be validated in this test. Given that a larger part of ALS heritability is due to rare gene mutations [62], we considered using a rare variant GWAS for validation, but the lack of publicly available summary statistics of well-powered studies of this kind limited this possibility.
We aimed to test the performance of our prediction method [24], previously validated for predicting drug side effects, in the context of gene discovery and prioritisation in ALS. However, we cannot claim that our method would perform better than alternative approaches as we did not perform any comparative assessment with respect to other techniques. The major difference between our method and most others is that ours is not based on a diffusion process such as random walks with restarts, and instead uses the local neighbourhood and an enrichment test to build the model. It will be important for future work to benchmark multiple methods across diseases to establish performance characteristics systematically across approaches, but that was not the aim of this work.
The genetics of ALS has proven challenging, and even though great progress has been made by generating and analysing large multi-omics datasets [36,[62][63][64][65][66][67][68][69][70][71], the causes of ALS in most patients (~85%) remain unexplained. Using machine learning models to leverage our current knowledge of ALS and other diseases could allow us to accelerate our progress in the understanding of the genetic causes of ALS and lead towards new avenues of treatment.

Research (NIHR) Biomedical Research Centres at South London & Maudsley and Guy s & St. Thomas NHS
Foundation Trusts, and part-funded by capital equipment grants from the Maudsley Charity (award 980) and Guy s & St. Thomas Charity (TR130505). The views expressed are those of the author(s) and not necessarily those of the NHS, the NIHR, King s College London, or the Department of Health and Social Care. We would like to thank people with MND and their families for their contribution to this project.

Conflicts of Interest:
The other authors declare no competing interests