# Accelerate the discovery of genetic variants in mitochondrial diseases with Variant prIOritization using Latent spAce

Justine Labory, Youssef Boulaimen, Jasmine Singh, Samira Ait-El-Mkadem Saadi, Véronique Paquis-Flucklinger, Sylvie Bannwarth, Silvia Bottini

<figure><img src="image_1.png" /></figure>

## To cite this version:

Justine Labory, Youssef Boulaimen, Jasmine Singh, Samira Ait-El-Mkadem Saadi, Véronique Paquis-Flucklinger, et al.. Accelerate the discovery of genetic variants in mitochondrial diseases with Variant prIOritization using Latent spAce. Briefings in Bioinformatics, 2025, 26 (6), pp.bbf612. 10.1093/bib/bbf612. hal-05469040

HAL Id: hal-05469040

https://hal.science/hal-05469040v1

Submitted on 21 Jan 2026

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Accelerate the discovery of genetic variants in mitochondrial diseases with Variant priOritization using Latent spAce

Justine Labory<fnref n="1" />,<fnref n="2" />, Youssef Boulaimen<fnref n="3" />, Jasmine Singh<fnref n="3" />, Samira Ait-El-Mkadem Saadi<fnref n="1" />,<fnref n="4" />, Véronique Paquis-Flucklinger<fnref n="1" />,<fnref n="4" />,

Sylvie Bannwarth<fnref n="1" />,<fnref n="4" />,†, Silvia Bottini ‡,<fnref n="2" />,†,*

## Abstract

Interpreting variants from whole-exome sequencing remains a major challenge, particularly for heterogeneous disorders such as mitochondrial diseases (MDs). To address this, we have developed Variant priLoritizationOn using Latent spAce (VIOLA), a pipeline designed to help find a diagnosis for complex cases. VIOLA uses a variational autoencoder to embed functional annotations into a low-dimensional space, followed by DBSCAN-based outlier detection to identify potential pathogenic variants. Filtering steps and phenotype integration via HPO terms are then applied. The VIOLA score (Vscore) combines variant outlierness, transcriptomic co-expression data, and MD-specific annotations. Two rankings are derived: the VIOLA rank (all variants) and the ARrank (variants compatible with autosomal recessive inheritance). The VIOLA Aggregated score (VAscore) merges Vscore with Exomiser's pathogenicity score. Applied to 20 patients (four diagnosed), VIOLA reduced the variant list by >99% and ranked causal variants within the top 5 using ARrank, outperforming existing methods. Overall, VIOLA is a patient-specific strategy for variant prioritization, helping to resolve challenging MD cases and uncover novel disease mechanisms.

**Keywords:** mitochondrial diseases; exome sequencing; variant prioritization; machine learning; multi-omics

## Introduction

Mitochondrial diseases (MDs) are rare and complex disorders caused by deficiencies in the mitochondrial respiratory chain, which is essential for cellular energy production through oxidative phosphorylation (OXPHOS) [1]. Primary mitochondrial diseases (PMDs) result from pathogenic variants in either mitochondrial DNA (mtDNA) or nuclear DNA (nDNA) directly affecting OXPHOS [2]. In contrast, secondary mitochondrial dysfunctions (SMDs) arise from defects in non-mitochondrial genes that indirectly impair mitochondrial functions [2-4]. This distinction is often difficult, as both can produce similar biochemical profiles [2].

Diagnosis of MDs is challenging due to genetic and clinical heterogeneity, with a wide range of symptoms and ages of onset [5]. Genetic transmission patterns add another layer of complexity: while mtDNA is exclusively maternally transmitted [6], nDNA follows Mendelian transmission [7], with around 80% of the known MD-associated genes being transmitted in an autosomal recessive (AR) manner [8].

Next Generation Sequencing (NGS) technologies such as Whole Exome Sequencing (WES) and Whole Genome Sequencing (WGS) greatly accelerated the discovery of disease-associated variants [9-13]. They identify tens of thousands of variants per patient, making it difficult to identify causal variants and requiring extensive curation by experts [11, 14]. To address this, variant prioritization tools have been developed [15], ranging from early annotation-based approaches such as CADD [16] and REVEL [17], to more advanced systems integrating phenotypic information such as LIRICAL [18], Exomiser [19-23], and Amelie [24]. At the same time, the ACMG/AMP guidelines [25] proposed by American College of Medical Genetics and the Association for Molecular Pathology provide standardized rules for the clinical interpretation of variants, and several tools automate their application [26-32]. Despite these advances, such strategies are often limited by incomplete data, such as lack of functional validation or segregation data, and prior ClinVar classification, leading to variants being classified as variants of uncertain significance (VUS). Importantly, none of the existing approaches are specifically tailored to MDs, particularly in complex cases that remain undiagnosed for

<footnote n="1">¹Université Côte d'Azur, Inserm U1081, CNRS UMR 7284, Institute for Research on Cancer and Aging Nice (IRCAN), 28, Ave de Valombrose, 06107 Nice Cedex 02, France</footnote>

<footnote n="2">²INRAE, Université Côte d'Azur CNRS, Institut Sophia Agrobiotech, 400 route des chappes, BP 167 06903 Sophia Antipolis Cedex, France</footnote>

<footnote n="3">³House of modeling, simulations and interactions, Université Côte d'Azur, Parc Valrose - 28, avenue Valrose Nice Cedex 2 06108, France</footnote>

<footnote n="4">⁴Department of Medical Genetics, Reference Centre for Mitochondrial Diseases Centre Hospitalier Universitaire (CHU) de Nice, 151 Route de Saint-Antoine, 06200 Nice, France</footnote>

<footnote>*Corresponding author. INRAE, Université Côte d'Azur, CNRS, Institut Sophia Agrobiotech, Sophia-Antipolis, France. E-mail: silvia.bottini@inrae.fr</footnote>

<footnote>‡Sylvie Bannwarth and Silvia Bottini authors are contributed equally.</footnote>
decades and involve unique variant combinations not captured by
standard pipelines.

Here, we introduce VIOLA (Variant PrIOritization using Latent
spAce), a novel multi-omics approach specific to MDs which
integrate genomic, transcriptomic and phenomic data and based
on the innovative hypothesis that since the disease-responsible
variant(s) are individual-specific and rare, therefore they should
have a unique combination of properties different from the rest
of the individual variants.

**Material and methods**

Patient cohort

We dispose of WES and RNA-sequencing data from 20 patients with suspicion of MDs, of which four with known responsible variant. The demographic details of the cohort are presented in Supplementary Table S1. For four patients, three causative genes were identified (the same gene for two patients since they are cousins), two of them were reported in Clinvar as pathogenic. The main characteristics of these variants can be found in Supplementary Table S2 and Section S1.

The processing steps for genomic and transcriptomic data are
detailed in Sections S2 to S5 of the Supplementary Material.

VIOLA workflow

VIOLA integrates genomics, phenomics, and transcriptomics to identify rare variants with peculiar combination of scores which could be potential candidates for diagnosis of very rare and particular MDs. The input of VIOLA is a VCF file upon annotation with Variant Effect Predictor (VEP) and CADD, and HPO terms linked to patient's phenotype. The output is a reduced ranked list of variants, which helps clinicians to speed up the identification of the causal variant during manual review.

The workflow is composed of three main blocks that are
described below (Fig. 1A).

Step A: variant features calculation

To characterize genetic variants, we used the CADD tool [16] which integrates multiple genomic features into one score to evaluate the deleteriousness of variants (Fig. 1A).

Therefore using CADD, we calculated 32 features to describe
variants: including six conservation scores from phastCons and
PhyloP tools [33], two pathogenicity scores from two functional
impact predictors: SIFT [34] and PolyPhen [35] and two scores
calculated by CADD [16]. We added eight epigenetic scores from
ENCODE catalog [36] and five splicing scores from MMSplice
tool [37]. Finally, the nine remaining scores were calculated
with BRAVO database [38] to include the genomic environment
of variants. All details about these features are described in
*Supplementary Table S3*.

Step B: techniques of dimensionality reduction

We used a variational autoencoder (VAE) to find complex nonlin-
ear combinations of variants features and reducing the dimen-
sionality.

The VAE was built with Tensorflow [39] and is composed of
one hidden layer in the encoder of dimension 25 and the latent
space which contains 16 latent features. Each layer of the encoder
or decoder is composed of the dense layer followed by a batch
normalization and the ELU activation function. To optimize the
model, a custom loss function was implemented, combining the
Mean Squared Logarithmic Error reconstruction loss and the Kull-
back-Leibler divergence to enforce latent space regularization.

We fed the VAE with the variants scores calculated in step A to
obtain 16 new features that are combinations of the 32 original
variants scores (Fig. 1A).

Step C: outlier detection

Since we reasoned that the responsible variant for MD suspicions
should have a very peculiar combination of parameters, this
should be seen as an outlier in the space of variants parameters.
We employed the Density-Based Spatial Clustering of Applications
with Noise (DBSCAN) algorithm as a first step toward the
identification of the outliers' variables [40]. DBSCAN is an unsupervised machine learning algorithm used for clustering spatial
data based on density. DBSCAN requires two parameters: the
maximum distance two points belonging to the same cluster
(epsilon, ε) and the minimum number of points required to form
a cluster (minPts).

For datasets with more than two dimensions, it is recom-
mended to set the minPts parameter to two times the number of
dimensions, therefore we settled to 32 [41].

To determine an appropriate ε, we used a data-driven approach based on the k-nearest neighbors (k-NN) method [42]. Specifically, for each point, we computed the average distance to its k-nearest neighbors, where k was set to $\sqrt{N}$ (with N being the number of variants). These distances were then sorted in ascending order, and the final ε value was chosen as the mean of the first *n* distances (where *n* was empirically set to 1500). Based on those parameters, DBSCAN divides the variants in two groups: the outliers and the non-outliers. We retrieve only the variants classified as outliers to continue the pipeline (Fig. 1A).

Step D: variant filtering

We applied a series of filters to variants retained from the previous
step (Fig. 1B). First, we removed all low-quality variants by exclud-
ing those with a quality score (QUAL field) below 50. The "QUAL"
field represents the Phred-scaled probability that a variant is truly
present at a given position. A score of 50 corresponds to an error
probability of 10⁻⁵, indicating very high confidence. This stringent
threshold was chosen to minimize false positives, as variants with
lower quality scores are often the result of sequencing errors or
low-confidence base calls [13, 43-45]. Next, since our study is
based on WES data, we focused exclusively on variants located
in protein-coding regions or which have an impact on protein-
coding genes. While WES primarily targets the exonic regions of
the genome where the majority of known pathogenic variants are
found, we also included intronic variants if they were predicted
to affect the function or regulation of a protein-coding gene. By
filtering for "protein-coding" variants, we reduced the data to
include only those variants that are likely to have functional
effects on gene products. Finally, we applied a final filter based
on variant allele frequency (VAF), which measures the proportion
of variant alleles at a genomic locus. For homozygous variants,
the VAF must be greater than 0.8, ensuring that the variant is
present in most sequencing reads at that position. Conversely, for
heterozygous variants, the VAF should fall between 0.3 and 0.7,
reflecting a balanced presence of both reference and alternative
alleles. Variants that do not meet these criteria are filtered out.

Step E: phenomics integration

To integrate the phenotype of the patient, we used the Human
Phenotype Ontology (HPO) [46]. HPO is a database that provides a
standardized vocabulary for describing phenotypic abnormalities
associated with human diseases.
<figure><img src="image_3.png" /><figcaption>Figure 1. Overview of the VIOLA pipeline. (A) The VIOLA pipeline begins with feature calculation using CADD from a VCF file. These 32 features (F1, F2...) are then transformed using a Variational Autoencoder (VAE) into 16 latent features (LF1, LF2...), followed by outlier detection with the DBSCAN algorithm. (B) The second stage of VIOLA focuses on filtering false positive variants. This is achieved by applying filters based on variant allele frequency (VAF), variant quality, and variant biotype. Next, the patient's phenotype is integrated using the Human Phenotype Ontology (HPO), retaining only variants that could potentially explain the observed phenotype. (C) In the final step, variants are ranked using different scoring methods. The first, the VIOLA score (Vscore), incorporates the Mahalanobis distance, transcriptomics data via a co-expression approach, and features related to mitochondrial diseases (MDs). The second, the VIOLA Aggregated Score (VAscore), combines the Vscore with the Exomiser score, a state-of-the-art tool for variant prioritization.</figcaption></figure>

To quantify the similarity between two phenotypic terms, we calculated the Jaccard Index for the sets of ancestors associated with the two HPO terms. We defined the Jaccard index of HPO terms as:

$$J(HPO_{patient}, HPO_{gene}) = \frac{|Ancestors_{HPO_{patient}} \cap Ancestors_{HPO_{gene}}|}{|Ancestors_{HPO_{patient}} \cup Ancestors_{HPO_{gene}}|}$$

The index gives values ranging from 0 to 1, where a value of 0 indicates that there are no shared elements between the sets (no similarity), and a value of 1 indicates that the sets are identical (complete similarity). Therefore, a high Jaccard index indicates greater similarity between the phenotypic features represented by the terms, while a low Jaccard index represents a small overlap (Fig. 1B).

### Step F: filtering variants based on the Mahalanobis distance

To identify variants that significantly deviate from the distribution of CADD scores, we computed the Mahalanobis distance [47], which measures the distance between a point and the mean of a multivariate distribution while considering correlations between variables. This approach allows for the detection of outliers based on multiple features simultaneously.

To do this, we first selected the 32 scores extracted in step A, then replaced missing values with the median of the corresponding score and finally computed Mahalanobis distance with the function *mahalanobis* from stats package [48].

To detect significantly distant variants, we calculated the P-value associated with the Mahalanobis distance using the $X^2$ distribution. Variants with $P < .001$ were considered outliers.

### Step G: transcriptomics data integration and definition of a first ranking score for variants: the VIOLA score

We defined a first score which includes multiple parameters to provide a first score to rank variants accordingly to the disease (Fig. 1C).

First, we included transcriptomics information leveraging MitoBook using a binary score. This score indicates whether a variant is associated with a gene belonging to a module enriched in mitochondrial genes, identified by co-expression analysis using WGCNA [49]. Variants linked to genes belonging to such modules are assigned a score of +1, while those outside these modules are assigned a score of -1. Details of this process are described in Supplementary Material in Section S6.

A second parameter used considers whether a variant is unique and specific to a patient, with the exception of cousins for whom we examined common variants not present in the rest of the
cohort. Variants meeting this criterion are assigned a score of +1, while non-unique variants are given a score of 0.

Prior knowledge regarding genes involved in MDs was also included. Referring to Stenton and Prokisch [8], 413 distinct genes are associated with MDs. Variants in these known genes are assigned a score of +1, as they are more likely to be implicated in the disease. Variants in genes not yet linked to MDs are assigned a score of 0.

Finally, although the Mahalanobis distance was already used to filter variants accordingly to the P-value associated to the measure, we also decided to include as additional term in the Vscore to consider the degree of significance provided by the statistical test.

To balance the contributions of the four components of the final Vscore, we assigned specific weights to each score. The Mahalanobis distance, transcriptomics, and uniqueness features were each assigned equal weights of 0.5, reflecting their comparable importance. In contrast, the score for known MD genes was assigned a lower weight of 0.01. This lower weighting accounts for the heterogeneity of MDs and the incompleteness of the current MD gene catalog.

The Vscore is calculated as the weighted sum of the Mahalanobis distance and the scores for transcriptomics, uniqueness, and known MD genes:

$$V \text{ score} = 0.5 \times (\text{Mahalanobis Distance} + \text{transcriptomics} + \text{uniqueness}) + 0.01 \times \text{known gene}.$$

Finally, the Vscore undergoes a normalization process using MinMax normalization to ensure that all scores fall within a standardized range, facilitating comparison across different variants. Once normalized, the variants are ranked based on their adjusted Vscore. Then, we apply a filtering step to retain only the most relevant variants by selecting those with a normalized Vscore greater than the third quartile (Q3), allowing us to focus on the top-ranking variants.

### Step H: definition of the autosomal recessive rank

Since most of the genes known to be involved in MDs are AR inherited, we created a rank specific for AR inheritance in which we incorporated the genotype of the patient called ARrank.

Variants are filtered according to their genotype. Only homozygous and compound heterozygous variants are compatible with AR mode of transmission, so we have kept all genes with homozygous variants and genes with at least two heterozygous variants.

The rank of the remaining variants is then recalculated still using the Vscore.

### Step I: definition of the VIOLA aggregated score

To create the VAscore, we calculated the arithmetic mean of the Vscore and the Exomiser score for each variant (Fig. 1C) (Details on how to run Exomiser are provided in Section S7 of the Supplementary Material). The Exomiser score ranges from 0 to 1 so each score contributes equally to the final VAscore as they have the same range.

### Step J: output of VIOLA

To improve the accessibility of our results, we have opted for a more interactive and user-friendly approach by providing the results in the form of an HTML file rather than as static tables, which can often be difficult to interpret. Interactive HTML tables bring together the most important information in a simplified format. These tables include key columns and scores, such as SIFT [34] and PolyPhen [35] values, MitoCarta database [50] and variant

frequencies from gnomAD [51] and 1000G [52]. This gives users rapid access to the most relevant data.

Databases such as dbSNP [53], gnomAD [51], OMIM [54], GeneCards [55], Genotype Tissue Expression (GTEx) [56] and Ensembl [57] are included to provide an overview of the gene or variant in question. By clicking on the links provided, users are directed to the relevant database entry in a new tab, offering immediate access to detailed information on the variant or gene in question.

## Implementation

VOLA is implemented in Python and R. The only two steps which are done in python are the use of the VAE and DBSCAN to find outliers. The rest of the pipeline is done in R.

## Results

### VIOLA: a new model to prioritize genetic variants

To develop a tool to perform variants prioritization for patients with suspected MD when the responsible variant is not present in common genes involved in MDs, we performed an extensive benchmark and selected the best performing tools. Our hypothesis is that since the disease-responsible variant is individual-specific and very rare, therefore it should have a unique combination of properties different from the other individual's variants. Hence, we assume that the putative disease variants for MD suspicions are outliers of each individual variants' distribution. VIOLA consists in three main blocks: feature calculation and feature extraction; variants filtering; variants prioritization. The main objective of VIOLA is to help reduce the list of variants and improve their ranking to be potential pathogenic compared to other variant prioritization tools to finally accelerate the identification of the responsible one(s).

We applied the pipeline on collections of rare variants from 20 patients with suspected MDs including 4 with known causal variants.

### Benchmark of feature calculation and feature extraction

To validate our approach and assess the robustness of our pipeline, we designed and tested several scenarios using data from the four patients for whom the causative variants are already known. These scenarios were designed to assess, firstly, the need for a dimensionality reduction technique and, secondly, which technique is the most effective (See more details in Section S7.1 of Supplementary Material). All the results are reported in Supplementary Table S4.

In the first scenario, we tested the pipeline without using the VAE (Supplementary Fig. S1A). Hence, we applied the DBSCAN clustering algorithm directly to the 32 scores calculated by CADD [16]. However, this approach proved ineffective, as all variants were assigned to the outlier group by DBSCAN. This result demonstrates that clustering the raw CADD scores without dimensionality reduction do not provide meaningful information for prioritizing potentially pathogenic variants in the context of MDs.

To further explore alternative approaches for dimensionality reduction, we set up two additional scenarios. Therefore, we replaced the VAE with two widely used dimensionality reduction methods: Principal Component Analysis (PCA) [58] (Supplementary Fig. S1B) and t-distributed Stochastic Neighbor Embedding (t-SNE) [59] (Supplementary Fig. S1C).

By feeding the PCA algorithm with the 32 scores, we extracted between 11 and 13 features for the four patients with known
<figure><img src="image_4.png" /></figure>

Figure 2. Proportion of retained variants at each step of the VIOLA pipeline across patients. The x-axis represents the 20 patients in the cohort, while the y-axis indicates the proportion of variants after each filter step of VIOLA relative to the initial set of variants. Each color corresponds to a different filtering step: outlier detection by DBSCAN (green), quality and allele frequency filters (orange), phenotype-based filtering using HPO terms (yellow), scoring with the VIOLA score (purple), and the final ranking with the VIOLA Aggregated score (pink). The stacked barplot illustrates how the number of candidate variants progressively decreases at each stage, highlighting the contribution of each filtering and scoring step in reducing the number of potential disease-causing variants.

causal variant, yielding between 11 633 and 90 958 variants classified as outliers by DBSCAN, indicating that PCA failed to adequately capture the critical features necessary for effective clustering (Supplementary Table S4). PCA's linear nature may not have been sufficient for modeling the complex relationships in the data.

Similarly, when t-SNE was used as feature extraction technique, the responsible variants were not classified as outliers. t-SNE failed to differentiate between the causative variants and the majority of benign variants.

Those results rule out PCA or t-SNE as alternative tools to VAE for feature extraction to be included in VIOLA.

The application of this step to our cohort filtered out 45% of variants (Fig. 2, Supplementary Table S5).

The integration of phenomics considerably reduces the number of variants classified as potential responsible for the disease.

After identifying the outlier group of variants using our pipeline, we applied a series of additional filters to refine the prioritized list of variants. These filters focused on three key criteria: the quality of the variants, the biotype of the genes affected by these variants and the VAF. Variants with low-quality scores were excluded to ensure the reliability of the analysis. Additionally, we restricted our analysis to variants located in protein-coding regions or other relevant biotypes known to play significant roles in disease mechanisms. The last filter concerned the VAF to remove variants that were inconsistent in terms of variant genotype and allele frequency. For more details on filters, see from Materials and methods. This first step resulted in a significant reduction in the number of variants to be screened, with an average of 57% of variants being eliminated at this stage (Fig. 2, Supplementary Table S5).

The next step focused on the identification of the variants relevant to the patient's phenotype. Therefore, we integrated phenotypic information from HPO database [46]. Clinicians have annotated the patients in the study with their correspondent HPO term(s) that exhaustively described the observed clinical phenotype. The challenge was to determine whether a given variant was linked to a gene associated with phenotypic features similar to those of the patient.

Therefore, we developed a novel approach based on the Jaccard index computation between the ancestors of HPO terms associated with the genes bearing variants and the HPO terms representing the patient's phenotype. This approach was iterated for all possible pairwise combinations of HPO terms associated with the variants and those of the patient. After calculating the

Jaccard index values, we arbitrarily established a threshold of 0.80 to retain only those variants associated with genes whose HPO terms showed significant similarity to the patient's phenotype. By applying this threshold, we ensured that the retained variants were likely to be functionally relevant to the clinical features observed in the patient. This is the most stringent step before the variants are ranked. On average, we excluded 85% of variants (Supplementary Table S5).

This score enables to include phenomics considerations into the selective process to identify the individuals' responsible variants.

Towards a strategy to include transcriptomics using co-expression network analysis.

To explore the relationships between genes and mitochondrial functions, we developed MitoBook a comprehensive catalog of gene co-expression patterns in patients with suspected MDs. Briefly, we applied a co-expression analysis to the transcriptomics data from the 20 patients from the cohort which we dispose. We applied WGCNA [49] to identify 12 distinct co-expression modules. These modules vary in size, ranging from 64 to 11 617 genes. Notably, the largest module, containing 11 617 genes, comprised almost half of the genes in the input expression matrix. Such a large module, lacking specificity, likely represents a background or generic expression pattern, it was excluded from further analyses to focus on modules with greater biological relevance. Similarly, the smallest module, consisting of only 64 genes, was also excluded.

Then, we investigated whether the 10 remaining modules were enriched in mitochondrial genes. To do this, we calculated a "propensity score" using the MitoCarta database [50] (see Supplementary Material, Section S6 for more details). This analysis enabled us to identify eight modules that show a significant enrichment in mitochondrial genes, ranging from 5.36 to 16.54%, indicating their propensity to be relevant in mitochondrial biology (Fig. 3A). The analysis was extended to additional nine databases containing information about the genes functions in mitochondrial functioning. The results are presented in the form of a heatmap showing, in rows, the different mitochondrial pathways we found, and, in columns, the modules found by WGCNA. The color scale represents the -log10 of the pvalue (Fig. 3B). The strongest enrichments were observed in "greenyellow" module which consists of 184 genes, 14.13% of which are mitochondrial genes. This module showed a significant overrepresentation of genes involved in mitochondrial ATP synthesis coupled to electron transport, as well as genes associated with mitochondrial
<figure><img src="image_5.png" /><figcaption>Figure 3. Mitochondrial gene proportion and pathway enrichment in co-expression modules. (A) The barplot shows the percentage of mitochondrial genes (according to the MitoCarta database) in each co-expression module. Modules highlighted in green represent those enriched in mitochondrial genes, while grey modules are not enriched. (B) The heatmap displays functional enrichment of the modules across five ontologies: Gene Ontology (GO) biological process (dark blue), GO cellular component, GO molecular function, human phenotype ontology (HPO), and Jensen diseases. The x-axis represents the co-expression modules, and the y-axis lists mitochondrial-related pathways identified. The color scale indicates the significance of enrichment (darker shades = stronger enrichment). The three modules with the highest enrichment in mitochondrial genes (A) also show strong enrichment for mitochondrial pathways (B), confirming their functional relevance.</figcaption></figure>

respiratory chain complex I and NADH dehydrogenase activity. This module groups together genes that play a crucial role in mitochondrial energy production and electron transport processes.

Overall, three modules with the highest propensity score from MitoCarta also achieved the highest score for the five other databases: GO Biological Process, GO Cellular Component, GO Molecular Function, HPO and Jensen Diseases. This convergence of evidence across multiple resources underlines the functional relevance of these three modules and highlights their potential importance in understanding mitochondrial processes and their links to patient phenotypes.

These results conducted to the definition of a score based on the guilty-by-association principle. Therefore, genes present in a module with high propensity score (>5.13), were assigned with a MitoBook score of 1 otherwise of -1.

The MitoBook score allows to consider an additional omics layer, the transcriptomic to be added to the overall estimation of the degree of pathogenicity.

## Evaluating the ranking efficacy of Vscore

To summarize all the calculated scores in one unique value to be used to prioritize potential pathogenic variants for MD suspicions, we first defined a novel comprehensive score, the VIOLA score (Vscore).

The Vscore incorporates multiple features, such as the Mahalanobis distance, the MitoBook score, and additional metrics related to MDs, to rank variants according to their likelihood of being pathogenic.

We chose Mahalanobis distance because it is a spatial measure and allows us to filter variants and detect outliers that

significantly differ from the distribution of selected features. We added Mitobook score to integrate transcriptomics data to better prioritize variants into genes that show similar expression profile to known mitochondrial genes. We also considered a score that reflects the uniqueness of the variant, since a variant that occurs more than once in our cohort (except for related patients) is less likely to be the responsible variant. Finally, we incorporated a score based on genes already known to be involved in MDs to add prior knowledge about the disease.

From the Vscore, we created two ranks, the Vrank which is general and the ARrank which includes genotype-specific filters. Exomiser score considers variant pathogenicity, frequency and the patient phenotype.

We developed two separate ranks because since most of the genes already involved in MDs have an AR inheritance, we wanted to give importance to genes compatible with AR inheritance by creating the ARrank while keeping a general rank with all the variants prioritized by VIOLA regardless of transmission mode (Vrank).

We evaluated the efficacy in ranking for both ranks of the Vscore on the four patients for whom the causal variants were already known, comparing the ranking obtained with Exomiser, the state-of-the-art tool to perform variant prioritization (for more details, see Section S7.3 of Supplementary Material). The results are summarized in Table 1.

For patients 3 and 17, the causative variant was ranked at the first place using both VIOLA ranks, it was ranked respectively at 9th and 17th positions for Exomiser (Fig. 4A and B). This indicates that the Vscore was able to identify the pathogenic variant as the most likely candidate, confirming its reliability in variant ranking.
<table><caption>Table 1. Ranks of responsible variants for diagnosed patients across different prioritization scores.</caption><thead><tr><td>Patient</td><td>Vrank</td><td>ARrank</td><td>VARank</td><td>Rank with Exomiser score</td></tr></thead><tbody><tr><td>3</td><td>1</td><td>1</td><td>1</td><td>9</td></tr><tr><td>4</td><td>7</td><td>3</td><td>3</td><td>18</td></tr><tr><td>9</td><td>5</td><td>4</td><td>1</td><td>19</td></tr><tr><td>17</td><td>1</td><td>1</td><td>1</td><td>17</td></tr></tbody></table>

The table presents the different ranking positions of the causal variants for the four diagnosed patients in the cohort. Columns indicate: the patient identifier (Patient), the rank of the causal variant obtained with the VIOLA score (Vrank), the rank when genotype is incorporated (ARrank), the rank obtained using the VIOLA Aggregated score (VARank) which combines VIOLA and Exomiser's scores. Lower rank values correspond to higher prioritization. This table shows that VIOLA outperforms Exomiser by ranking causal variants in the top 5 positions, whereas Exomiser ranks them in the top 20.

A

Patient 3

<figure><img src="image_6.png" /></figure>

B

Patient 17

<figure><img src="image_7.png" /></figure>

C

Patient 4

<figure><img src="image_8.png" /></figure>

D

Patient 9

<figure><img src="image_9.png" /></figure>

Figure 4. Comparison of responsible variant ranks using different scoring methods for diagnosed patients. Radar plots compare the ranking performance of different prioritization strategies for the responsible variants of four diagnosed patients: patient 3 (A), patient 17 (B), patient 4 (C), and patient 9 (D). Each axis represents a ranking method, VIOLA rank (Vrank), Autosomal Recessive rank (ARrank), VIOLA Aggregated rank (VARank), and Exomiser rank, and we plotted the top 20 positions only. The distance from the center indicates the position of the responsible variant (lower values represent better ranks). These plots highlight differences in ranking efficiency across methods, showing that VIOLA consistently places the causal variants among the top candidates, outperforming Exomiser.

For patient 4, the Vscore reported the responsible variant in the 7th position for the Vrank while it was 18th with Exomiser. However, if the ARrank is used, the rank improved significantly, placing the variant in the third position (Fig. 4C).

Similarly, for patient 9, the responsible variant achieved a rank of fifth when using the Vrank. With the ARrank, the rank improved further, elevating the variant to the fourth position while with Exomiser, the variant was ranked 19th (Fig. 4D).

Therefore, we defined an integrated score considering both Exomiser score and Vscore into a single score, the VIOLA Aggregated Score (VAScore) (Fig. 4, Table 1). By exploiting the unique strengths of each tool, the VAscore offers a more robust and

accurate ranking system, improving the identification of pathogenic variants in complex rare disease cases.

The results highlight the added value of integrating genotype-specific data to refine the prioritization process, particularly in cases where the responsible variant may not stand out as the top candidate based solely on other features.

## Validation by conducting shuffle analysis

To validate the VIOLA pipeline, since we do not dispose of another cohort, we performed a shuffled analysis on the feature values associated to each variant as the input of VIOLA. We designed an
approach to randomize the feature values to assess the impor-
tance of these scores for the ranking of the causal variants.
Using data from the four patients with known causal variants,
we shuffled the features values 100 times while preserving their
original distribution. After each shuffle, we applied the entire
VIOLA pipeline (see Section S7.2 of Supplementary Material for
more details).

At the DBSCAN clustering level, an interesting pattern emerged: the clustering of variants in the input files remained consistent with the results obtained using the original CADD scores. However, a significant effect was observed when calculating Vscore and Vrank. In 28%-69% of the shuffled analyses, depending on the patient, the causal variant was filtered out and therefore not classified by VIOLA. For the remaining analyses where the causal variant was retained, we examined its rank distribution to assess the impact of shuffling, and all the results are summarized in Supplementary Table S6.

For patients 3 and 17, shuffling seems to have a less significant
effect on causal variant ranking than expected. The Vrank ranged
from 1 to 12 for patient 3 and from 1 to 10 for patient 17. The
ranking obtained for the responsible variants from the original
variants parameters was preserved in only four out of 100 of the
shuffled analyses for both patients.

For patients 4 and 9 the shuffling impacted dramatically the ranking. For patient 9, the causal variant reached a maximum Vrank of 54 and a median of 19 and the original rank was only reproduced twice out of 72, while the highest Vrank of causal variant of patient 4 was 1184.

These results underline the crucial importance of variants
parameters and their combination in our pipeline. The shuffled
analyses reinforced the ability of VIOLA to accurately prioritize
and identify the causal variant. This underlines that the combina-
tion of parameters from original scores are fundamental for the
stability and accuracy of the rankings.

Results on the other patients of the cohort

We extended the application of VIOLA to the undiagnosed
patients of our cohort. In absence of the knowledge about
the responsible variants, we compared the performances of
VIOLA and Exomiser regarding the prioritization of variants
in mitochondrial-related genes. Therefore, we compared the
proportion of variants located in genes listed in the MitoCarta
database for both tools. Across the 20 patients in the cohort,
VIOLA consistently identified a higher proportion of variants
related to mitochondrial genes than Exomiser (Fig. 5), suggesting
that VIOLA is more specifically tailored to prioritize variants that
are functionally relevant to mitochondrial biology. By focusing on
features and pathways associated with mitochondrial function,
VIOLA can better capture disease-relevant variants that may be
overlooked by more general tools such as Exomiser.

Discussion

In this study, we presented VIOLA, a novel variant prioritization
tool that integrate machine learning and multi-omics data to
facilitate the identification of causal variants in the challenging
context of MDs. By combining genomics, transcriptomics, and
phenomics in a single framework, VIOLA addresses several lim-
itations of existing prioritization methods.

The diagnosis of MDs is particularly complex due to their
heterogeneity, the difficulty of distinguishing primary between
SMDs, and the incomplete annotation of mitochondrial genes.
As a result, some patients remain in diagnostic stalemate for
decades. In some cases, pathogenic variants lie in genes not yet

linked to mitochondrial function, making them especially difficult
to detect.

Recent benchmarking studies report that the most effective
variant prioritization tools are those integrating both genomic
and phenotypic information [15, 60–63]. For example, Kelly *et al.*
identified 28 tools that accept a VCF file and HPO terms as input
[15]. However, many of these tools are no longer maintained, are
not freely accessible, or lack local and programmatic access, an
essential feature for processing sensitive VCF data [15].

Overall, benchmarks reported Exomiser, LIRICAL, AMELIE, and xRare as the top-performing tools. However, each of these has limitations: AMELIE does not support local execution [15], xRare is not compatible with GRCh38 data [60], and Exomiser, while considered the most mature by bio.tools [63], was not designed to capture the complex patient-specific variant signatures observed in MDs. Recent approaches, such as ACMG/AMP-guideline-based probabilistic scoring [31] or ensemble frameworks like EWE3 [64], which combines Exomiser, xRare, and DeepPVP [65] are emerging, but they remain limited to supervised strategies relying on predefined knowledge.

VIOLA differs from these approaches by introducing an
unsupervised machine learning framework specifically tailored
to MDs, integrating not only genomic and phenotypic data but
also transcriptomic co-expression and MD-specific features into
a unified prioritization score.

The aim of VIOLA is to reduce the list of variants for manual
review and improving ranking performance compared to existing
tools. In our data, VIOLA ranked the causal variants of all four
diagnosed patients within the top five positions using the ARrank.
Importantly, three of these variants were located in genes not cur-
rently annotated as mitochondrial, underlining VIOLA's ability to
identify novel disease mechanisms in the context of MDs beyond
known mitochondrial gene sets such as MitoCarta [50] or the list
of genes already known to be involved in MDs [8]. This highlights
the tool's potential not only for accurate variant prioritization but
also as an exploratory resource to uncover previously unknown
factors in MDs. VIOLA can therefore help diagnose complex or
atypical cases, where conventional approaches may prove insuf-
ficient. In addition, VIOLA provides three complementary scores
that can be integrated with pathogenicity scores from existing
tools or with ACMG/AMP guidelines based on the patient's profile.

However, several limitations remain. Our evaluation was
limited by the small size of the cohort, as genomic data from
patients with suspected or diagnosed MDs are both rare and
strictly regulated. This limitation partly motivated the design
of VIOLA as a patient-specific tool. Furthermore, our aim was to
develop a tool to be employed in diagnostic context where usually
small cohorts are mainly available. Despite this constraint,
VIOLA succeeded in improving variant prioritization compared
to Exomiser and even ranked the same causal variant in first
position in two related patients, without prior knowledge of their
relationship.

Another limitation is that our approach is based on WES data
and therefore excludes noncoding variants, which may play key
roles in gene regulation or splicing. While WGS could address this
[66], the large volume of generated data considerably increases
the complexity of the analysis, making it difficult to identify the
responsible variant among millions, of potential candidates.

Future work should focus on extending the scope of the
analysis to WGS data and incorporating functional validation
experiments to confirm the pathogenicity of the highest-ranked
variants. In addition, the weighting scheme for different Vscore
features could be further optimized, possibly using larger and
more diverse datasets.
<figure><img src="image_10.png" /></figure>

Figure 5. Proportion of candidate variants in mitochondrial genes prioritized by VIOLA and Exomiser. Stacked barplots illustrating the proportion of candidate variants located in genes listed in the MitoCarta database across 20 individual patients. For each patient, two bars are shown, corresponding to the variants prioritized by VIOLA (upper panel) and Exomiser (lower panel). The y-axis indicates the percentage of variants assigned to mitochondrial (green) or non-mitochondrial (grey) genes. The x-axis represents individual patients of the studied cohort. This figure shows that VIOLA consistently identified a higher proportion of variants in mitochondrial genes than Exomiser.

## Conclusion

In conclusion, VIOLA represents an exploratory tool to help diagnose complex cases of patients with suspected MDs. By combining advanced computational methods with prior knowledge and different omics layer, VIOLA improves the identification of causal variants, paving the way for more effective diagnostics and personalized treatments for patients with suspected MDs and potentially other rare diseases. Of note, VIOLA can easily be extended to other rare diseases with minor adjustments to the Vscore, thus broadening its scope of application. Further development and validation of VIOLA on larger cohorts and diversified datasets will be essential to establish its generalizability and robustness in different clinical and research settings.

### Key Points

* VIOLA is a new tool designed to help identify pathogenic variants on the assumption that pathogenic variants exhibit a unique combination of properties in each individual.

* VIOLA combines unsupervised learning via a variational autoencoder, outlier detection using DBSCAN, and phenotype-based filtering to effectively reduce the number of candidate variants.

* VIOLA integrates genomic, phenomic and transcriptional co-expression data and mitochondrial-specific annotations into a unified score (Vscore), with an additional ARrank adapted to autosomal recessive inheritance models.

* VIOLA shows robust performance with real data. Applied to a cohort of patients with suspected mitochondrial disease, VIOLA systematically ranked the causal variants among the top 5 candidates with the ARrank.

## Acknowledgements

The authors are grateful to the Genotoul bioinformatics platform Toulouse Midi for providing computing and storage resources.

## Author contributions

Justine Labory (Data Curation, Conceptualization, Formal analysis, Methodology, Software, Writing—original draft, Writing—review & editing), Youssef Boulaimen (Methodology, Software), Jasmine Singh (Methodology, Software), Samira Ait-El-Mkadem Saadi (Writing—original draft, Writing—review & editing), Veronique Paquis (Writing—review & editing), Sylvie Bannwarth (Supervision, Writing—review & editing), and Silvia Bottini (Conceptualization, Supervision, Writing—original draft, Writing—review & editing)

## Supplementary data

Supplementary data is available at *Briefings in Bioinformatics* online.

Conflict of interest: None declared.

## Funding

This work was supported by the Agence Nationale pour la Recherche (ANR 21-PMRB-0006) and the French Foundation For Rare Disease (Fondation Maladies Rares). This work was supported by the French government, through the UCA JEDI Investments in the Future project managed by the National Research Agency (ANR) under reference number ANR-15-IDEX-01.

## Data availability

VIOLA is implemented in R and python and is available in the following GitHub repository: https://github.com/Plant-Net/VIOLA.

## References

1. Gorman GS, Chinnery PF, DiMauro S. et al. *Mitochondrial diseases Nature Reviews Disease Primers* 2016;2:1–22.

2. Parikh S, Karaa A, Goldstein A. et al. Diagnosis of ‘possible’ mitochondrial disease: An existential crisis. *J Med Genet* 2019;56: 123–30. https://doi.org/10.1136/jmedgenet-2018-105800.
3. Niyazov DM, Kahler SG, Frye RE. Primary mitochondrial disease and secondary mitochondrial dysfunction: Importance of distinction for diagnosis and treatment. *Molecular Syndromology* 2016;7:122-37. https://doi.org/10.1159/000446586.

4. Mancuso M, McFarland R, Klopstock T. et al. International workshop:: Outcome measures and clinical trial readiness in primary mitochondrial myopathies in children and adults. Consensus recommendations. *Neuromuscul Disord* 2017;27:1126-37. https://doi.org/10.1016/j.nmd.2017.08.006.

5. Munnich A, Rustin P. Clinical spectrum and diagnosis of mitochondrial disorders. *Am J Med Genet* 2001;**106**:4-17. https://doi.org/10.1002/ajmg.1391.

6. et al. Maternal inheritance of mammalian mitochondrial DNA. *Nature* 1974;251:536-8. https://doi.org/10.1038/251536a0.

7. Alston CL, Rocha MC, Lax NZ. et al. The genetics and pathology of mitochondrial disease. *J Pathol* 2017;241:236-50. https://doi.org/10.1002/path.4809.

8. Stenton SL, Prokisch H. *Genetics of Mitochondrial Diseases*. Identifying mutations to help diagnosis. *EBioMedicine*, Elsevier, 2020. 102784, **56**. https://doi.org/10.1016/j.ebiom.2020.102784.

9. Wortmann SB, Koolen DA, Smeitink JA. et al. Whole exome sequencing of suspected mitochondrial patients in clinical practice. *J Inherit Metab Dis* 2015;**38**:437-43. https://doi.org/10.1007/s10545-015-9823-y.

10. Souche E, Beltran S, Brosens E. et al. Recommendations for whole genome sequencing in diagnostics for rare diseases. *Eur J Hum Genet* 2022;**30**:1017-21. https://doi.org/10.1038/s41431-022-01113-x.

11. Davis RL, Kumar KR, Puttick C. et al. Use of whole-genome sequencing for mitochondrial disease diagnosis. *Neurology* 2022; **99**:e730-42. https://doi.org/10.1212/WNL.000000000200745.

12. Dewey FE, Grove ME, Pan C. et al. Clinical interpretation and implications of whole-genome sequencing. *JAMA* 2014;**311**: 1035-45. https://doi.org/10.1001/jama.2014.1717.

13. Riley LG, Cowley MJ, Gayevskiy V. et al. The diagnostic utility of genome sequencing in a pediatric cohort with suspected mitochondrial disease. *Genet Med* 2020;**22**:1254-61. https://doi.org/10.1038/s41436-020-0793-6.

14. Caspar S, Dubacher N, Kopps A. et al. Clinical sequencing: From raw data to diagnosis with lifetime value. *Clin Genet* 2018;**93**: 508-19. https://doi.org/10.1111/cge.13190.

15. Kelly C, Szabo A, Pontikos N. et al. Phenotype-aware prioritisation of rare mendelian disease variants. *Trends Genet* 2022;**38**: 1271-83. https://doi.org/10.1016/j.tig.2022.07.002.

16. Kircher M, Witten DM, Jain P. et al. A general framework for estimating the relative pathogenicity of human genetic variants. *Nat Genet* 2014;**46**:310-5. https://doi.org/10.1038/ng.2892.

17. Ioannidis NM, Rothstein JH, Pejaver V. et al. REVEL: An ensemble method for predicting the pathogenicity of rare missense variants. *Am J Hum Genet* 2016;**99**:877-85. https://doi.org/10.1016/j.ajhg.2016.08.016.

18. Robinson PN, Ravanmehr V, Jacobsen JO. et al. Interpretable clinical genomics with a likelihood ratio paradigm. *Am J Hum Genet* 2020;**107**:403-17. https://doi.org/10.1016/j.ajhg.2020.06.021.

19. Robinson PN, Köhler S, Oellrich A. et al. Improved exome prioritization of disease genes through cross-species phenotype comparison. *Genome Res* 2014;**24**:340-8. https://doi.org/10.1101/gr.160325.113.

20. Smedley D, Jacobsen JOB, Jäger M. et al. Next-generation diagnostics and disease-gene discovery with the exomiser. *Nat Protoc* 2015;**10**:2004-15. https://doi.org/10.1038/nprot.2015.124.

21. Bone WP, Washington NL, Buske OJ. et al. Computational evaluation of exome sequence data using human and model organism

phenotypes improves diagnostic efficiency. *Genet Med* 2016;**18**: 608-17. https://doi.org/10.1038/gim.2015.137.

22. Cipriani V, Pontikos N, Arno G. et al. An improved phenotype-driven tool for rare mendelian variant prioritization: Benchmarking exomiser on real patient whole-exome data. *Gen* 2020;**11**:460. https://doi.org/10.3390/genes11040460.

23. Vestito L, Jacobsen JO, Walker S. et al. Efficient reinterpretation of rare disease cases using exomiser. *NPJ Genom Med* 2024;**9**:65.

24. Birgmeier J, Haeussler M, Deisseroth CA. et al. AMELIE speeds mendelian diagnosis by matching patient phenotype and genotype to primary literature. *Sci Transl Med* 2020;**12**:eaaau9113.

25. Richards S, Aziz N, Bale S. et al. Standards and guidelines for the interpretation of sequence variants: A joint consensus recommendation of the American College of Medical Genetics and Genomics and the Association for Molecular Pathology. *Genet Med* 2015;**17**:405-23. https://doi.org/10.1038/gim.2015.30.

26. Tavtigian SV, Greenblatt MS, Harrison SM. et al. Modeling the ACMG/AMP variant classification guidelines as a Bayesian classification framework. *Genet Med* 2018;**20**:1054-60. https://doi.org/10.1038/gim.2017.210.

27. Nicora G, Limongelli I, Gambelli P. et al. CardioVAI: An automatic implementation of ACMG-AMP variant interpretation guidelines in the diagnosis of cardiovascular diseases. *Hum Mutat* 2018;**39**: 1835-46. https://doi.org/10.1002/humu.23665.

28. Xavier A, Scott RJ, Talseth-Palmer BA. TAPES: A tool for assessment and prioritisation in exome studies. *PLoS Comput Biol* 2019;**15**:e1007453. https://doi.org/10.1371/journal.pcbi.1007453.

29. Scott AD, Huang K-L, Weerasinghe A. et al. CharGer: Clinical characterization of germline variants. *Bioinformatics* 2019;**35**: 865-7. https://doi.org/10.1093/bioinformatics/bty649.

30. Li Q, Zhao K, Bustamante CD. et al. Xrare: A machine learning method jointly modeling phenotypes and genetic evidence for rare disease diagnosis. *Genet Med* 2019;**21**:2126-34. https://doi.org/10.1038/s41436-019-0439-8.

31. Nicora G, Zucca S, Limongelli I. et al. A machine learning approach based on ACMG/AMP guidelines for genomic variant classification and prioritization. *Sci Rep* 2022;**12**:2517. https://doi.org/10.1038/s41598-022-06547-3.

32. Wang Z, Zhao G, Zhu Z. et al. VarCards2: An integrated genetic and clinical database for ACMG-AMP variant-interpretation guidelines in the human whole genome. *Nucleic Acids Res* 2024;**52**:D1478-89. https://doi.org/10.1093/nar/gkad1061.

33. Siepel A, Bejerano G, Pedersen JS. et al. Evolutionarily conserved elements in vertebrate, insect, worm, and yeast genomes. *Genome Res* 2005;**15**:1034-50. https://doi.org/10.1101/gr.3715005.

34. Ng PC, Henikoff S. Predicting deleterious amino acid substitutions. *Genome Res* 2001;**11**:863-74. https://doi.org/10.1101/gr.176601.

35. Adzhubei IA, Schmidt S, Peshkin L. et al. A method and server for predicting damaging missense mutations. *Nat Methods* 2010;**7**: 248-9. https://doi.org/10.1038/nmeth0410-248.

36. Abascal F, Acosta R, Addleman NJ. et al. Expanded encyclopaedias of DNA elements in the human and mouse genomes. *Nature* 2020;**583**:699-710. https://doi.org/10.1038/s41586-020-2493-4.

37. Cheng J, Nguyen TYD, Cygan KJ. et al. MMSplice: Modular modeling improves the predictions of genetic variant effects on splicing. *Genome Biol* 2019;**20**:48. https://doi.org/10.1186/s13059-019-1653-z.

38. Taliun D, Harris DN, Kessler MD. et al. Sequencing of 53,831 diverse genomes from the NHLBI TOPMed program. *Nature* 2021;**590**:290-9. https://doi.org/10.1038/s41586-021-03205-y.

39. Abadi M, Agarwal A, Barham P. et al. TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems. OSDI'16: Proceedings
of the 12th USENIX conference on Operating Systems Design and Implementation. 2015.

40. Ester M, Kriegel H-P, Sander J. et al. A density-based algorithm for discovering clusters in large spatial databases with noise 1996;96: 226-31.

41. Sander J, Ester M, Kriegel H-P. et al. Density-based clustering in spatial databases: The algorithm gdbscan and its applications. Data Min Knowl Disc 1998;2:169-94. https://doi.org/10.1023/A:1009745219419.

42. Rahmah N, Sitanggang IS. Determination of Optimal Epsilon (Eps) Value on Dbscan Algorithm to Clustering Data on Peatland Hotspots in Sumatra, Vol. 31. IOP Conference Series: Earth and Environmental Science, Volume 31, Workshop and International Seminar on Science of Complex Natural Systems 9-10 October 2015, Bogor, Indonesia, 2016, 012012.

43. Kremer LS, Bader DM, Mertes C. et al. Genetic diagnosis of mendelian disorders via RNA sequencing. Nat Commun 2017;8:15824. https://doi.org/10.1038/ncomms15824.

44. Krumm N, Turner TN, Baker C. et al. Excess of rare, inherited truncating mutations in autism. Nat Genet 2015;47:582-8. https://doi.org/10.1038/ng.3303.

45. Cooke DP, Wedge DC, Lunter G. Benchmarking small-variant genotyping in polyploids. Genome Res 2022;32:403-8. https://doi.org/10.1101/gr.275579.121.

46. Gargano MA, Matentzoglu N, Coleman B. et al. The human phenotype ontology in 2024: Phenotypes around the world. Nucleic Acids Res 2024;52:D1333-46. https://doi.org/10.1093/nar/gkad1005.

47. De Maesschalck R, Jouan-Rimbaud D, Massart DL. The Mahalanobis distance. Chemom Intell Lab Syst 2000;50:1-18. https://doi.org/10.1016/S0169-7439(99)00047-7.

48. R Core Team R. R: A language and environment for statistical computing 2013.

49. Langfelder P, Horvath S. WGCNA: An R package for weighted correlation network analysis. BMC Bioinformatics 2008;9:1-13. https://doi.org/10.1186/1471-2105-9-559.

50. Rath S, Sharma R, Gupta R. et al. MitoCarta3.0: An updated mitochondrial proteome now with sub-organelle localization and pathway annotations. Nucleic Acids Res 2021;49:D1541-7. https://doi.org/10.1093/nar/gkaa1011.

51. Chen S, Francioli LC, Goodrich JK. et al. A genomic mutational constraint map using variation in 76,156 human genomes. Nature 2024;625:92-100. https://doi.org/10.1038/s41586-023-06045-0.

52. 1000 Genomes Project Consortium. A global reference for human genetic variation. Nature 2015;526:68. https://doi.org/10.1038/nature15393.

53. Sherry ST, Ward M, Sirotkin K. dbSNP—Database for single nucleotide polymorphisms and other classes of minor genetic variation. Genome Res 1999;9:677-9. https://doi.org/10.1101/gr.9.8.677.

54. Hamosh A, Scott AF, Amberger JS. et al. Online mendelian inheritance in man (OMIM), a knowledgebase of human genes and genetic disorders. Nucleic Acids Res 2005;33:D514-7. https://doi.org/10.1093/nar/gki033.

55. Stelzer G, Rosen N, Plaschkes I. et al. The GeneCards suite: From gene data mining to disease genome sequence analyses. Curr Protoc Bioinformatics 2016;54:1-30.

56. GTEx Consortium, Ardlie KG, Deluca DS. et al. The genotype-tissue expression (GTEx) pilot analysis: Multitissue gene regulation in humans. Science 2015;348:648-60.

57. Dyer SC, Austine-Orimoloye O, Azov AG. et al. Ensembl 2025. Nucleic Acids Res 2025;53:D948-57. https://doi.org/10.1093/nar/gkae1071.

58. Wold S, Esbensen K, Geladi P. Principal component analysis. Chemom Intell Lab Syst 1987;2:37-52. https://doi.org/10.1016/0169-7439(87)80084-9.

59. Van der Maaten L, Hinton G. Visualizing data using t-SNE. J Mach Learn Res 2008;9:2579-2605.

60. Jacobsen JO, Kelly C, Cipriani V. et al. Phenotype-driven approaches to enhance variant prioritization and diagnosis of rare disease. Hum Mutat 2022;43:1071-81. https://doi.org/10.1002/humu.24380.

61. Yuan X, Wang J, Dai B. et al. Evaluation of phenotype-driven gene prioritization methods for mendelian diseases. Brief Bioinform 2022;23:bbac019. https://doi.org/10.1093/bib/ bbac019.

62. Tosco-Herrera E, Muñoz-Barrera A, Jaspez D. et al. Evaluation of a whole-exome sequencing pipeline and benchmarking of causal germline variant prioritizers. Hum Mutat 2022;43:2010-20. https://doi.org/10.1002/humu.24459.

63. Licata L, Via A, Turina P. et al. Resources and tools for rare disease variant interpretation. Front Mol Biosci 2023;10:1169109.

64. Fan Y, Zhou Y, Liu H. et al. Improving variant prioritization in exome analysis by entropy-weighted ensemble of multiple tools. Clin Genet 2023;103:190-9. https://doi.org/10.1111/cge.14257.

65. Boudellioua I, Kulmanov M, Schofield PN. et al. DeepPVP: Phenotype-based prioritization of causative variants using deep learning. BMC Bioinformatics 2019;20:65.

66. Mattick JS, Dinger M, Schonrock N. et al. Whole genome sequencing provides better diagnostic yield and future value than whole exome sequencing: The integration of genome sequencing with clinical records and data from the internet of things will transform health care. Med J Aust 2018;209:197-9. https://doi.org/10.5694/mja17.01176.