Generating artificial dilution series

Motivation and requirements.

We assume that we have a set of genes, each with a varying number of GO annotations. This set of annotated genes will be referred to as the truth set, T. T is expected to be a table that contains gene names in the first column and the associated GO term identifiers in the second column. We assume the existence of an ontology, such as GO. The ontology links smaller, more precise child terms to more broadly defined ancestral terms. Here, if a term t_child, is a child term for t_par, then all genes that belong to t_child will also belong to t_par (i.e. t_child ⊆ t_par). These child–ancestor connections form a directed acyclic graph, where every GO term, with the exception of the root term, has a varying number of parent terms [13]. Genes can be mapped to GO terms at varying levels in the ontology, depending on how detailed our level of knowledge is.

Artificial Dilution Series creates several data sets, each representing functional classifications for genes in T. These data sets, referred to here as Artificial Prediction Sets (or AP sets), imitate a GO prediction set from a GO classifier. They are comprised of a list of triples, (gene, GO, sc), where sc is the artificial prediction score. Each AP set is created artificially, with a defined proportion of correct and erroneous annotations. AP set creation does not require a real classifier, nor the actual gene sequences. Furthermore, each AP set can also be defined as a table, where rows represents different annotation entities and columns represent genes, GO terms and classifier scores.

The workflow for the ADS pipeline is shown in Fig 1 and outlined in pseudocode in S1 Text. The first step is to create two GO prediction sets, the Positive Prediction set, P_pos, and the Negative Prediction set, P_neg. These steps use a signal model to perturb the ground truth GO annotations and a noise model to generate erroneous annotation predictions. Note that both P_pos and P_neg have a table structure, similar to the AP set.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 1. ADS workflow.

The ADS pipeline iterates through specified signal levels (e.g. from 100% to 0%) at each level creating a collection of artificial prediction sets (AP sets). Each AP set is a union of a set containing only false annotations (negative set) and a set containing a controlled fraction of true annotations (positive set). AP sets are compared against the correct annotations using an Evaluation Metric (EvM). The results can be shown graphically by plotting EvM scores for AP sets at each signal level as boxplots. These reveal how stable an EvM is against random variation introduced at different signal levels and whether it can measure the amount of signal retained in different AP sets. Finally, EvM performance is quantified with rank correlation.

https://doi.org/10.1371/journal.pcbi.1007419.g001

Step 1: Create the positive set.

The Positive set, P_pos, is created from a copy of T in two stages: First, we draw a random integer, N_shift, between 0 and the size of P_pos, and select N_shift random rows from P_pos. Next, the GO term, t_corr, of every selected row is replaced with a semantically similar term. This term is selected from the k nearest parents of t_corr. We refer to this procedure as shifting to semantic neighbours and it forms the basis of our signal model. Figure 1 in supplementary text S2 Text demonstrates this step. We have tested k = (2, 3, 4) and found results to be similar.

Second, we introduce a percentage of errors, ϵ, by switching GO terms between genes for a controlled fraction of annotation lines (see Fig 2 for a toy example). We refer to this as the permutation of the positive set. Here a random pair of annotations, a = (gene_a, t_a) and b = (gene_b, t_b), where gene_a ≠ gene_b and t_a ≠ t_b, is selected from P_pos. If t_a is sufficiently far away in the GO structure from every GO term annotated to gene_b and vice versa for t_b and the annotations of gene_a, then the GO labels are exchanged. The distance between t_x and t_y is defined as the Jaccard correlation between the two sets created from the ancestral classes of t_x and t_y (see S1 Text). In our experiments we used a threshold, th_noise, of 0.2 and required that the correlation was lower than th_noise. This step represents our noise model.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 2. ADS permutation example.

In this toy example we start with a set of correct annotations using only two GO terms A and B, taken from truth set T. Initially, in the upper part the Artificial Prediction Set (APS) matches perfectly to T. This would give us 100% correct predictions. Next, the permutation step, shown in the middle, switches GO terms for a pair of genes. This increases false positives and false negatives by 1 and decreases true positives and true negatives by 1. We obtain as a result an APS with a known level of errors (0.002 in this example). The permutation is repeated with other gene pairs and other GO terms until the required noise level is obtained.

https://doi.org/10.1371/journal.pcbi.1007419.g002

The permutation process is repeated until the required percentage of errors, ϵ, have been introduced. This adds a controlled fraction of false positives and false negatives among the predictions. We refer to this fraction as the ADS noise level and the remaining fraction of correct positive annotations as the ADS signal level (1—ADS noise level).

Motivation for signal and noise models

We have two reasons for allowing t_signal in the signal model to move to nearby terms in the ontology. First, real-life GO prediction tools often predict nearby terms as these tend to be correlated with t_corr. Second, omitting these changes would favour metrics that focus on t_corr in the evaluation and ignore the ontology structure during evaluation, which would bias our analysis results.

The noise model, embedded in the Positive Set, aims to represent the worst type of error a classifier could make as a) t_noise is clearly separated from t_corr in the ontology and b) prediction scores for erroneous predictions, sc_noise, come from the same distribution as correct predictions. This means that a threshold on the sc scores is unable to separate signal and noise predictions, ensuring that noise cannot be excluded from the evaluation. An EvM should consider our noise predictions as negative cases and score the data set accordingly. Furthermore, generating noise by permuting the GO classes from T keeps many other variables constant, such as average class size and average distance to the root node of predicted terms. The comparison of evaluation metrics is therefore not biased towards metrics that consider unrelated variables.

The noise that is added with lower prediction scores is motivated by trends in classifier results, where lower sc corresponds to lower prediction quality. Note that most metrics take this into consideration by using varying thresholds on the sc values. If this noise was excluded, there would be no benefit in using a threshold on AP set results as an EvM could simply average over all prediction results.

Testing evaluation metrics with ADS

EvMs are tested by scoring the AP sets and comparing the results to the truth set, T. This gives us an output matrix: (1) where Metric is a given evaluation metric, columns, j = [1, 2, ‥k], correspond to different ADS signal levels, rows i = [1, 2‥N] represents the different repetitions with the same ADS signal level and APS_ij is the i^th AP set obtained at signal level j. A good metric should have low variance within one ADS signal level (within Output columns) and a clear separation between the different ADS signal levels (between Output columns). We test this in two ways: A) visually by plotting Output columns as separate boxplots at different ADS noise levels and B) numerically by calculating the rank correlation, RC, between Output matrix and ADS signal levels. Boxplots give an overview of the performance of the tested metric. RC monitors how closely EvM is able to predict the correct ranking of the ADS signal levels. One could also consider the linear correlation, but that would favour linear correlation over the correct rank.

False positive data sets

In addition to AP sets, we also generated several prediction sets where all genes were annotated with the same set of non-informative GO terms or with randomly chosen GO terms. We wanted to be able to identify evaluation metrics that give high scores to prediction sets containing no signal. We refer to these sets as False Positive Sets (FP sets). We included the following FP sets:

1. Each gene has the same large set of the smallest GO terms. We call this the all positive set.
2. Each gene has the same small set of the largest GO terms, where sc is taken from the term frequency in the database. We call this the naïve set.
3. Each gene has a set of randomly selected GO terms. We call this the random set.

The motivation for the all positive set is to test the sensitivity of the evaluation metric to false positive predictions (low precision). Note that a set where each gene has all GO terms would be impractical to analyse. Therefore, we decided to select a large set of the smallest GO terms, each having many ancestral GO terms. The random set shows the performance when both precision and recall are very low and should therefore always have a low score. The naïve set tests if the evaluation metric can be misled by bias in the base rates to consider the results meaningful. The naïve set uses the same list of GO classes, ranked in decreasing order of size, for every analysed gene. The naïve set was originally proposed by Clark and Radivojak [21] and it has been used extensively in CAFA competitions [2]. They refer to the naïve set as naïve method or naïve predictor. Note that all these FP sets are decoupled from the annotated sequences and contain no real information about the corresponding genes.

There is a parameter that needs to be considered for FP set generation: the number of GO terms reported per gene. We selected this to be quite high (800 GO terms). This emphasises the potential inability of the metric to separate the noise signal in the FP sets from the true signal in predictions output by classifiers. We used three false positive sets:

1. Naïve-800: all genes assigned to the 800 most frequent GO terms. Here larger GO classes always have a better score.
2. Small-800: all genes assigned to the 800 least frequent GO terms. Here smaller GO classes always have a better score.
3. Random-800: all genes assigned by a random sample of 800 GO terms (sampling without replacement).

GO frequencies for FP set sampling was estimated using the UniProt GO annotation (ftp://ftp.ebi.ac.uk/pub/databases/GO/goa/UniProt/, 01.2019).

In the visual analysis, we plot the score for each FP set as a horizontal line over the boxplot of AP set scores. A well-performing evaluation metric is expected to score all FP sets close to AP sets with minimal or no signal (signal ∼ 0), i.e. as a horizontal line at the lowest region of the plot (see Fig 3).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 3. Combining FP sets with ADS results.

Here F_max results from False Positive (FP) sets are combined with F_max boxplots from ADS. FP results are added as horizontal lines to the visualisation. A robust metric would score FP sets near the lowest boxplot. An FP signal estimate is obtained by comparing each FP result with medians of each signal level. This is shown in figure with vertical arrows. Resulting signals are here limited to between 0 and 1.

https://doi.org/10.1371/journal.pcbi.1007419.g003

In our numerical analysis, we summarise FP sets by 1) Finding the approximate point where the horizontal FP set line crosses the curve created by the medians of each signal level (see Fig 3), 2) choosing the corresponding signal level from the x-axis as the FP score for the FP set in question and 3) selecting the maximum over all FP sets as the final FP Score (referred to later as FPS). Note that in step 2 we limit the signal between zero and one. In step 3 we focus on the worst result over all tested FP sets. A good EvM should not be sensitive towards any of the FP sets, making the worst signal (= highest FP score) the best measure for the insensitivity towards these biases.

Evaluation metrics

We selected representatives from various types of EvMs commonly used in AFP evaluation (see supplementary text S2 Text). We also included modifications of some EvMs. Note that there is a large variety of EvMs available, especially in the machine learning literature [8, 10], and we cannot cover them all. The tested EvMs can be grouped into three families: rank-based metrics, GO semantic similarity-based metrics and group-based metrics.

We evaluated two rank-based metrics: area under ROC curve (AUC-ROC) and area under precision-recall curve (AUC-PR). Both strategies capture the amount of false positive and false negative errors at different threshold levels. We evaluated three semantic similarity-based metrics: Lin [22], Resnik [23] and Ancestor Jaccard (AJacc). AJacc is a novel method taking the Jaccard correlation between two sets of GO terms (i.e. the ancestor terms for both GO terms). Group-based methods compare a set of predicted GO terms to the ground truth. From this group, we evaluated SimGIC [20], Jaccard correlation (Jacc) and S_min [24]. We also evaluated F_max, one of the most popular metrics in machine learning.

Most of the metrics compared here, have been used actively in AFP evaluation (see Table 1 in supplementary text S2 Text). Although AUC-PR has not been used directly, it is related to the popular precision-recall plots. All metrics except rank-based metrics can be calculated using different thresholds over prediction scores. We selected the maximum score as the final score for each metric over all possible threshold values.

The semantic similarity metrics we tested return a matrix, M_g, of similarity values for each analysed gene, g. Here each m_i,j in M is a similarity score between the predicted GO term, i, and the correct GO term, j, and is calculated using the selected semantic similarity metric. Altogether M_g compares the set of correct GO terms for g, S_corr, and the set of predicted GO terms for g, S_pred, with each other. Here S_corr is a fixed set, whereas the S_pred is altered, including GO predictions for g that have stronger prediction scores than the acceptance threshold, th. We iterate through thresholds from the largest to the smallest values of the prediction scores resulting in a series of S_pred sets for each gene with the smallest set corresponding to the largest th and the largest to the smallest th. Next, we generate M_g for every S_pred and summarise M_g into a single score with a procedure described below. From threshold-wise scores we then select the best (for most metrics the maximum) score.

When semantic similarity is used as an evaluation metric, M_g needs to be summarised into a single score for each gene at each tested threshold. As there is no clear recommendations for this, we tested six alternatives:

1. A. Mean of matrix. This is the overall similarity between all classes in S_corr and S_pred
2. B. Mean of column maxima. This is the average of best hits in M_g for classes in S_corr.
3. C. Mean of row maxima. This is the average of best hits in M_g for classes in S_pred.
4. D. Mean of B and C.
5. E. Minimum of B and C.
6. F. Mean of concatenated row and column maxima.

Methods A and D have been used previously [19]. A is considered to be weaker than D because A is affected by the heterogeneity of GO classes in S_corr (see Supplementary Text S2 Text and [19]). Methods B and C represent intentionally flawed methods, used here as negative controls. B is weak at monitoring false positive predictions and C is weak at monitoring false negatives. D aims to correct B and C by averaging them, but is still sensitive to outliers. Therefore, we propose methods E and F as improvements to method D. E monitors the weaker of B and C, and is therefore more stringent at reporting good similarity. F combines two vectors used for B and C. This emphasises the longer of two vectors as the size of S_pred is altered. Section 5.7 in our supplementary text S2 Text motivates these summation methods in more detail. Finally, the obtained scores for each gene at each th threshold are averaged and the maximum over all th values is reported.

Most EvMs can be further modified by using the same core function with different Data Structuring. By data structuring we refer to the gene centric, term centric or unstructured evaluations. In Gene Centric evaluation (GC), we evaluate the predicted set of GO terms against the true set separately for each gene and then summarise these values with the mean over all genes. In Term Centric evaluation (TC) we compare the set of genes assigned to a given GO term in predictions against the true set separately for each GO term and take the mean over all GO terms. In UnStructured evaluation (US) we compare predictions as a single set of gene–GO pairs, disregarding any grouping by shared genes or GO terms. In total we tested 37 metrics (summarised in S1 Table). Further discussion and metric definitions are given in S2 Text.

Data sets and parameters

All metrics were tested on three annotation data sets: the evaluation set used in CAFA 1 [2], the evaluation set from the MouseFunc competition [4] and a random sample of 1000 well-annotated genes from the UniProt database (ftp://ftp.ebi.ac.uk/pub/databases/GO/goa/UniProt/, 01.2019). These data sets represent three alternatives for evaluating an AFP tool (see Introduction). We used annotations from the beginning of 2019 with all data sets. It should be noted that the annotation densities of the tested data sets vary over time, potentially affecting the ranking of metrics.

These data sets have different features, shown in S2 Table. MouseFunc is the largest with 45,035 direct gene–GO annotations over 1634 genes. After propagation to ancestors, this expands to 291,627 GO annotations, composed of 13,691 unique GO terms. This gives a matrix with 1.42% of cells containing annotations. The CAFA data set had the smallest number of genes, with 7,092 direct gene–GO annotations over only 698 genes. After propagation, there are 61,470 GO annotations, mapping to 6,414 unique GO terms. This gives a matrix with 1.1% of cells containing annotations. The Uniprot data had the smallest number annotations, with 5,138 GO annotations over 1000 genes. These gave 46,884 annotations to 5,138 classes after propagation with 0.9% of cells in the resulting matrix containing annotations. All three data sets are intentionally small, as this is a frequent challenge in AFP evaluation. Furthermore, the number of annotations is skewed across both genes and GO terms in all data sets (see S2 Table).

We tested all metrics with three values of k (semantic neighbourhood size): k = 2, 3, 4. RC and FPS scores for these runs are listed in S4 Table (k = 2), S5 Table (k = 4) and S6 Table (k = 3). We also tested th_noise = 0.8 with k = 3. Metric scores for all AP and FP sets are available in S1, S2 and S3 Files.

Stability of ADS results

ADS data generation includes two main parameters: the size of the accepted neighbourhood, k, while shifting to semantic neighbours and a threshold for accepting noise classes, th_jacc. Here we show that results are stable as these values vary. We first compared the results with k = 2, 3, 4. We set th_jacc = 0.2. We created ADS and FP results with each value of k and calculated the Pearson and rank correlations of the results for k = 2 vs. 3 and 2 vs. 4. As these were already highly correlated, it was unnecessary to proceed any further. Analysis was done separately on RC and FP results and repeated on each tested data set.

The results, shown with sub-table A in S3 Table, show that variations appear insignificant. The Pearson correlation ranges between 0.9944 and 0.9999. Similarly, rank correlation ranges between 0.9935 and 1. This suggests that results are stable against small variations in k. We selected results with k = 3 for all further analyses.

We also tested the effect of Jaccard threshold, th_jacc, in our noise model. We compared our standard value th_jacc = 0.2 to th_jacc = 0.8. We performed a similar evaluation as with values of k (see sub-table B in S3 Table). The minimum and maximum Pearson correlation was 0.9991 and 0.9999, respectively, and the minimum and maximum rank correlation was 0.9945 and 0.9997, suggesting that our results are stable against th_jacc

ADS identifies differences between evaluation metrics

We compared the following evaluation metrics: US AUC-ROC, F_max, S_min (referred as S_min2 in figure), Resnik with summary method D and A, and Lin with summary method D to see if ADS can find performance differences. Visualisations for these metrics are shown in Fig 4 (for exact values see S6 Table).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 4. Six popular evaluation metrics compared with ADS.

Visual analysis of ADS results for six evaluation metrics: US AUC-ROC, F_max, S_min2, Resnik score A and D, and Lin score D, obtained with CAFA data. Scores for AP sets at each ADS signal level are shown as boxplots and scores for FP sets as horizontal lines. RC value shows rank correlation of the AP sets with ADS signals and FPS shows the highest signal from FP sets. RC should be high and FPS should be low. Note the drastic differences between methods. We discuss these in the main text. We flipped the sign of S_min results for consistency. ADS signal is plotted on the x-axis and the evaluation metric, shown in headings, is plotted on the y-axis. Horizontal lines for FP sets are: Naive = red line, All Positive = blue line and Random = Green line.

https://doi.org/10.1371/journal.pcbi.1007419.g004

Results with CAFA data show drastic differences between evaluation metrics. We see a clear separation across ADS signal levels in the boxplots for F_max (RC = 0.982) and S_min (RC = 0.985), and next best separation is for US AUC-ROC. Three of the semantic similarity measures performed poorly: they have unstable distributions in the obtained RC scores. Resnik in particular, is shown to have the worst separation (RC = 0.808 and 0.811 for summary A and D, respectively). Our results show that it is impossible to infer the original signal levels using these particular semantic similarity measures.

It could be argued that comparing random data and data with a strong signal should be sufficient for metric evaluation. Our results, however, show that in most plots the first (signal = 1) and last (signal = 0) boxplots do not fully reveal how good the evaluation metric is. For Lin score and Resnik score, for example, two extreme boxplots might be acceptable, assuming clear separation. However, the intermediate levels reveal considerable variability.

With FP sets visualised with horizontal lines, we see different methods failing. Now AUC-ROC has the worst performance as it ranks the FP sets as equally good as signal = 1. Note, that FP sets do not convey any real annotation information, but represent biases related to the GO structure. Altogether, only Resnik with method A has excellent performance with the FP signal and only S_min has good performance in both tests. These results demonstrate that neither ADS nor FP sets are sufficient to expose all weaknesses in isolation. They clearly show orthogonal views on the evaluated methods, both of which are important.

ADS confirms known flaws in AUC-ROC

Here we focus on the performance of AUC methods, using area under receiver operating characteristic curve (AUC-ROC) and area under precision-recall curve (AUC-PR) with the CAFA data set. AUC-ROC has been used extensively in AFP evaluation, whereas AUC-PR has not. We look at the unstructured (US), gene centric (GC) and term centric (TC) versions of both AUC methods.

AUC-ROC has been criticised for being a noisy evaluation metric with sensitivity to sample size and class imbalance [8, 25]. ADS demonstrates that AUC-ROC metrics correlate reasonably well with signal level, but fail with FP sets (Fig 5 and S6 Table). These results are expected, as FP sets monitor sensitivity to class imbalance. Surprisingly, two of our of three AUC-ROC methods rank the FP sets as good as signal = 1 (Fig 5, top row). This means that near-perfect prediction results could not be distinguished from false positives using this metric.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 5. AUC metrics compared using ADS.

We compared standard AUC (AUC-ROC) and AUC for precision-recall curve (AUC-PR) with CAFA dataset. We test unstructured (US), gene centric (GC) and term centric (TC) versions. US AUC-ROC and GC AUC-ROC exhibited exceptionally poor performance with the FP sets. TC AUC-PR showed good performance with both ADS and FP sets. RC and FPS are explained in previous figure text. ADS signal is plotted on the x-axis and the evaluation metric, shown in headings, is plotted on the y-axis.

https://doi.org/10.1371/journal.pcbi.1007419.g005

The poor FP set performance of some AUC metrics is caused by the identical treatment of small and large GO classes in US and GC analysis. However, a positive result for a large class is more probable from a random data set. This is worsened by the extreme class size differences in GO. This problem was corrected in the MouseFunc competition [4] by dividing GO classes to subsets based on class sizes and in CAFA 2 [3] by analysing each GO term separately, as in TC AUC-ROC. Indeed, the TC analysis displayed in the last column performs well with the FP sets. This is replicated across all data sets.

When AUC-ROC and AUC-PR are compared, while we see that FP sets have a weaker signal with US AUC-PR and GC AUC-PR, the problem persists. In the correlation scores, however, AUC-PR metrics consistently outperform AUC-ROC metrics. The difference between AUC-ROC and AUC-PR might be explained by the definition of AUC-ROC (see supplementary S2 Text). Since TC AUC-PR assigns low scores to FP sets and achieves high correlation, it is our recommended method with TC AUC-ROC in second place. Furthermore, these results argue strongly against the usage of US and GC versions of AUC metrics.

ADS reveals drastic differences between semantic summation methods

We compared the following GO semantic similarities: Resnik, Lin and Ancestor Jaccard (AJacc), each coupled with six different semantic summation methods, A-F (see Methods). The results are shown in Fig 6 and S6 Table. Fig 6 is organised so that the columns contain the compared semantic similarities and rows contain the summation methods. Ideally, we want a combination that delivers good separation in boxplots and ranks FP sets next to AP sets with signal = 0.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 6. Semantic similarities with different summation methods.

We present every combination of three semantic similarity-based methods (Resnik, Lin, AJacc) in columns and six semantic summation methods (A-F, see Methods) in rows. We show again the results for CAFA dataset. The summation methods have a bigger impact on performance than the actual metric. The novel summation methods, E and F, outperform the previous standards, A and D. ADS signal is plotted on the x-axis and the evaluation metric combined with a given summation method is plotted on the y-axis.

https://doi.org/10.1371/journal.pcbi.1007419.g006

This comparison shows the power of ADS. Our first observation is that the choice of summation method has a bigger impact than which semantic measure is used. Methods A and C fail at separating signal levels in ADS and methods B and D fail at ranking FP sets next to the signal = 0. The only methods that perform well in both ADS and FP tests are E and F. Furthermore, the intentionally flawed methods, B and C, were found to be weak by either ADS or FP tests.

Our second observation is that Lin similarities exhibit marginally better performance than the other metrics tested, as shown in RC values. While these results might suggest that Lin similarity summarised with methods E or F performs the best, this is dependent on data set (see below).

S_min and SimGIC have similar good performance

The third group of metrics we tested were methods that score sets of predicted GO terms. Some methods, like SimGIC and SimUI [19], calculate the Jaccard correlation between predictions and the ground truth, whereas S_min calculates a distance between them [26]. Most of these metrics include Information Content (IC) weights [19] in their function. We have two alternatives for IC, but limit the results here to IC weights proposed in the S_min article [26] (details on two versions of information content, ic and ic2, are discussed in the S2 Text). We still show results for both ic and ic2 in later analyses. We further include gene centric and unstructured versions of the following functions:

1. SimGIC. Gene centric original version (weighted Jaccard)
2. SimGIC2. Unstructured modified version of SimGIC
3. S_min1. S_min original version
4. S_min2. Gene centric version of S_min
5. GC Jacc. Gene centric unweighted Jaccard (original SimUI)
6. US Jacc. Unstructured version of SimUI (unstructured Jaccard)

The results are shown in Fig 7 and S6 Table. Both S_min and SimGIC show very high correlation with the ADS signal level. However, in FP tests, SimGIC has better performance. Also, US Jacc shows quite good performance, but GC Jacc is weaker. Introducing gene centric terms in S_min2 resulted in a notable increase in FP scores in all three data sets. SimGIC2, our unstructured version of SimGIC, showed a minor improvement in correlation with all three data sets. Overall, the unstructured versions (S_min1, SimGIC2 and US Jacc) outperformed the gene centric versions (S_min2, SimGIC, GC Jacc) in correlation values, FP tests or both. Performance of S_min and SimGIC metrics across different data sets was surprisingly similar. Our recommendation is to use S_min1 or SimGIC2. Alternatively, one can use US Jacc, despite weaker performance than S_min and SimGIC2, when a simpler function is required.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 7. S_min and SimGIC metrics.

We compared two versions of S_min, SimGIC and Jaccard each with Uniprot data. Performance was similar with all metrics, with the exception of GC Jacc which had a weaker correlation. S_min1 and SimGIC2 are slightly better methods. ADS signal is plotted on the x-axis and the evaluation metric, shown in headings, is plotted on the y-axis. We flipped the sign of the S_min metrics for consistency.

https://doi.org/10.1371/journal.pcbi.1007419.g007

Metric performance varies between data sets

We wanted to compare all the discussed evaluation metrics to see which of them performed best overall. We plotted RC against FPS for all metrics tested on the CAFA data set (see Fig 8). High quality metrics appear in the upper left corner of the plot. Metrics that fail with FP tests appear in the right side of the plot. Metrics failing in ADS correlation will be lower in the plot.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 8. Top metrics for CAFA data set.

Metrics are shown with different colours. F_max is black, AUC-ROC and AUC-PR are in red, Jacc, S_min and SimGIC are in green and semantic similarity metrics are in blue. Colours are varied between different metrics in the same group. The figure uses the following abbreviations: AJ = AJacc, F = F_max, G = SimGIC, J = Jacc, L = Lin, PR = AUC-PR, R = Resnik, ROC = AUC-ROC and Sm = S_min. Abbreviations are combined with summation or data structuring method. For example, Resnik F is R-F, Lin A is L-A, TC AUCROC is TC-ROC, ic2 SimGIC is ic2-G and ic S_min is ic-S. High quality metrics appear in the upper left corner (high correlation and low False Positive Signal). This clustered part is zoomed in the lower panel for clarity.

https://doi.org/10.1371/journal.pcbi.1007419.g008

The upper panel in Fig 8 shows big differences in performance. The semantic similarity metrics (shown in blue) show drastic scatter on the y-axis and AUC metrics (red) show drastic scatter on the x-axis, in agreement with results from earlier sections. Clearly, usage of unstructured and gene centric Area Under Curve methods (GC AUC-ROC, US AUC-ROC, US AUC-PR, GC AUC-PR) and semantic similarities based on signal summation methods A, B, C and D cannot be recommended. The metrics with the best performance (lower panel in Fig 8) are Lin E, Lin F, TC AUCPR, two versions of SimGIC2 and S_min2.

When we repeated the same analysis with the other two data sets, the rankings are different. We show only the best-scoring metrics in Fig 9 for MouseFunc and Uniprot. The performance of all metrics is shown in S1 and S2 Figs. In MouseFunc, the best performing group is the SimGIC—S_min methods with S_min2 as the best method. Lin E and Lin F are also among the top performers. TC AUC-PR, one of previous best performers, has weaker performance. S_min improved its performance, while TC AUC-PR got weaker. The Uniprot data set had differing results as well, with TC AUC-PR, Resnik E, Resnik F and TC AUC-ROC showing the best performance. The biggest drop among the top performing metrics is with Lin semantics (Lin E and Lin F).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 9. Top metrics for MouseFunc (top panel) and Uniprot (lower panel) data sets.

Here we repeat the lower plot from Fig 8 on the other two data sets. The performance of individual EvMs varies between data sets. Colouring and labels are explained in Fig 7.

https://doi.org/10.1371/journal.pcbi.1007419.g009

Semantic similarities had surprising behaviour in these experiments. Lin similarity (with methods E and F) was clearly the best performing semantic similarity-based metric in CAFA and MouseFunc, with Resnik as the weakest. However, on Uniprot data, Resnik improves its performance, making it the best semantic similarity method, while Lin’s performance drops significantly. We are unable to explain these differences.

Another surprising finding is that we see higher FP signals for TC AUC-PR than expected, particularly in the MouseFunc data. As term centric analysis evaluates each GO class separately, it should be insensitive to biases created by favouring larger or smaller GO classes. We will explain this in the next section.

Note that these data sets differed from each other in size and density, with Uniprot being the smallest and sparsest data set, CAFA being intermediate and MouseFunc data being the largest and densest data set. We see the following tendencies in the results: a) Lin semantic performs well until the size and sparsity of the Uniprot data is reached, b) term centric AUC-ROC and AUC-PR perform well except with denser and larger mouse data and c) simGIC methods appear to perform consistently. Elucidating the exact reasons for these performance differences would require further studies. Altogether, this suggests that selecting the best evaluation metric is data-dependent, showing the necessity of ADS-style tests.

Selecting overall best evaluation metrics

Results above showed variation in the ordering of metrics across the datasets. Still, it is important to see if we can select for AFP method comparisons metrics with best overall performance. Therefore, we collected the results of top performing metrics for each data set into Table 1, showing RC (Rank Correlation) and FPS (False Positive Set) values. In addition, we show the same results for popular methods with weaker performance. Again, a good metric should have high RC values and low FP values. Therefore, we show the five highest scores in RC tests in each column with bold and five lowest scores in underlined italics. In FP tests, we similarly show five lowest scores in bold and five highest scores in underlined italics. Excellent metrics should have bold numbers and no italics.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 1. Summary of results for best performing and widely-used metrics.

Here we show RC (Rank Correlation) and FP (False Positive) results for the best performing methods. We also show same results for some widely-used metrics. Good metrics should have a high RC score and low FP scores. Rec column shows our selected recommendations (See text for details). The five best results in each column are shown in bold. The five weakest results in each column are shown with underlined italics. Metrics that fail a given test are highlighted in red (see text for details). Note how methods in lower block show consistent weak performance either in RC or FP tests.

https://doi.org/10.1371/journal.pcbi.1007419.t001

Furthermore, we wanted to highlight the consistent high performance of metrics by setting a threshold on RC, t_RC, and threshold on FPS, t_FPS, and requiring that RC > t_RC and FPS < t_FPS. We highlight the failures of these tests in Table 1 in red. For RC values, we observed that the top performing metrics would be selected with threshold t_RC = [0.94, 0.98] and therefore set t_RC = 0.95. For FP values, we found that the top performing metrics would be selected with t_FPS = [0.10, 0.25] and set t_FPS = 0.16. There is a trade-off when setting these thresholds, as different values can emphasise either RC or FPS. However, both tests should be monitored for metric selection. In conclusion, a good metric should have as little red shading as possible and high performance in both RC and FP tests. Our supplementary table S6 Table shows these results for all the tested metrics. We selected the best performing metrics from supplementary table to the upper portion of Table 1.

Conclusion

Science needs numerous predictive methods, making the comparison and evaluation of these methods an important challenge. The central component of these comparisons is the Evaluation Metric. In the case of Automated Function Predictors the selection of the evaluation metric is a challenging task, potentially favouring poorly performing methods. In this article we have described the Artificial Dilution Series (ADS) which is, to our knowledge, the first method for testing classifier evaluation metrics using real data sets and controlled levels of embedded signal. ADS allows for simple testing of evaluation metrics to see how easily they can separate different signal levels. Furthermore, we test with different False Positive Sets representing data sets that the evaluation metric might mistakenly consider meaningful. Our results show clear performance differences between evaluation metrics. We also tested several modifications to existing evaluation metrics, some of which improved their performance. This work provides a platform for selecting evaluation metrics for specific Gene Ontology data sets and will simplify the development of novel evaluation metrics. The same principles could be applied to other problem domains, inside and outside of the biosciences.