Wei2GO: weighted sequence similarity-based protein function prediction

Wei2GO input preparation

Wei2GO requires two input files:

1. Diamond or BLASTP hits against the UniProtKB database in tabular format (--outfmt 6).

2. HMMScan hits against the Pfam-A database as a table output (--tblout).

These files can either be produced separately, or be incorporated with Wei2GO via its Snakemake pipeline (Mölder et al., 2021).

Wei2GO algorithm

Each protein matched with a UniProtKB reference via BLAST or DIAMOND is annotated with GO terms associated to these reference proteins in the GOA database. Consecutively, each protein matched with a Pfam domain via HMMScan is annotated with GO terms associated to this domain in the pfam2go database.

For every protein v, a GO term i is given a weight W(v,i) by taking the sum of the -log10 e-value E of all UniProtKB and Pfam sequence similarity hits j associated with this GO term, where e-values above 1 are filtered out:

(1) $$W\left(v,i\right)=\sum _{j\in Hits\left(v,i\right)}-log10\left(E_j\right),$$

where Hits(v,i) is the set of BLAST and HMMScan sequence aligment hits to v annotated with GO term i, and $E_j$ is the e-value of the sequence alignment j. Unique to Wei2GO compared to Argot is that terms annotated to a UniprotKB protein with an evidence code other than IEA (Inferred from Electronical Annotation) are given a weight multiplier:

(2) $$W\left(v,i\right)=\sum _{j\in Hits\left(v,i\right)}-log10\left(E_j\right)\cdot log2\left(P\right),$$

where P is the ratio between the number of GO terms annotated with an IEA vs GO terms annotated with an other evidence code in the GOA database. Given that the amount of GO terms annotated with an IEA code is substantially higher than those with non-IEA codes, P is larger than 1 and results in higher weights when based on non-IEA instead of IEA evidence codes.

In accordance with the GO Directed Acyclic Graph (DAG) structure all weights associated to child terms of a GO term are added to its weight, with the rationale that if a protein is annotated with a certain GO term it by default includes all its parent terms. This has the effect of assigning more confidence in GO terms higher up the DAG, i.e. annotations that are more generic and less prone to errors.

The total weight associated with each protein-GO term pair is used to calculate the Internal Confidence (InC) score by taking the relative weight of the GO term compared to the weight of the root node (3):

(3) $$\mathrm{I}\mathrm{n}\mathrm{C}\left(\mathrm{v},\mathrm{i}\right)={\displaystyle \frac{\mathrm{W}\left(v,i\right)}{W\left(root\right)}.}$$

The InC score is used to calculate a Group Score (GS) for each GO term i associated to a protein, by taking the sum of the InC of all GO terms j annotated to a protein v, if GO term j has a semantic similarity of 0.7 or higher with GO term i:

(4) $$GS\left(v,i\right)={\sum }_{j\epsilon annotations\left(v\right)}InC\left(v,j\right)\cdot I\left(Sim\left(i,j\right)\ge 0.7\right),$$

where I(.) is the indicator function which is equal to 1 if its argument is true and 0 otherwise, and annotations(v) is the set of GO annotations for protein v. Note that each GO term can have a unique GS score, as opposed to Argot where one singular GS score is defined for a group of similar GO terms. The value of Sim (semantic similarity) is calculated using Lin’s formula (Lin, 1998) as:

(5) $$Sim\left(i,j\right)={\displaystyle \frac{2\cdot \mathrm{I}\mathrm{C}\left(k\right)}{IC\left(i\right)+IC\left(j\right)},}$$

where k is the first joint parent term between GO term i and j, and the IC (Information Content) is calculated as:

(6) $IC\left(i\right)=-\mathrm{l}\mathrm{o}\mathrm{g}10(P\left(i\right),$

where P(i) is the fraction of proteins annotated to term i. Finally, a total score (TS) for each GO term annotated to a protein v is calculated as:

(7) $$TS\left(v,i\right)=IC\left(i\right)\cdot {\displaystyle \frac{\mathrm{I}\mathrm{n}\mathrm{C}\left(v,i\right)}{GS\left(v,i\right)}\cdot \mathrm{W}\left(v,i\right).}$$

The TS score is similar to that of Argot2 and Argot2.5, which combines the scores calculated by the equations above. However, Argot introduces an additional non-cumulative InC score, i.e. Eq. (3) without adding the weights of child terms, to remove bias of scores towards terms which have a higher probability of being predicted, i.e. a lower IC. This is done in an attempt to reduce bias towards predicting GO terms that are generic at the cost of GO terms that are highly specific, as reflected by their IC. Wei2GO instead operates under the assumption that GO terms which are easier to predict due to any reason, including having a lower IC, should have this reflected in their score. Any possible negative effects of this bias towards terms with lower IC affecting terms with higher IC will be reflected in the ‘S_min’ evaluation metric, which measures performance proportional to these IC’s (See ‘Test set evaluation’ section).

Test set generation

To assess the performance of Wei2GO it is compared to Argot2.5 and its predecessor Argot2, with the addition of DeepGOPlus (Kulmanov & Hoehndorf, 2020) which acts as a reference point. All proteins in the GOA database created between 01-01-2018 and 01-06-2020 were extracted. Only GO terms with the evidence code EXP (Inferred from Experiment), IDA (Inferred from Direct Assay), IPI (Inferred from Physical Interaction), IMP (Inferred from Mutant Phenotype), IGI (Inferred from Genetic Interaction), IEP (Inferred from Expression Pattern), TAS (Traceable Author Statement), and IC (Inferred by Curator), were retained. All molecular function terms ‘GO:0005515’ (protein binding) were removed from the annotation set if it was the only molecular function term annotated to the protein, to remove a large amount of bias due to many proteins having this generic term as their only molecular function term. This same approach was used in the CAFA3 test set generation.

The test set was expanded by including all terms according to their GO DAG structure using the ‘is_a’ and ‘part_of’ relationships for parent terms. The final test set consists of 4,067 proteins with 97,054 assigned GO terms of which 1,922 proteins have 59,697 biological process terms, 651 proteins have 4,785 molecular function terms, and 2,387 proteins have 32,572 cellular component terms.

Test set evaluation

Precision and recall is used for the assessment of Wei2GO. For a threshold t, the precision pr and recall rc are calculated per protein v, $T_v$ is the set of true GO terms for a protein v, $P_v$ is the set of predicted GO terms for a protein v, m(t) is the number of proteins with at least one GO term predicted with a score above the threshold, and n is the total number of proteins:

(8) $$p{r}_v\left(t\right)={\displaystyle \frac{\left|P_v\cap T_v\right|}{\left|P_v\right|},}$$

(9) $$r{c}_v\left(t\right)={\displaystyle \frac{\left|P_v\cap T_v\right|}{\left|T_v\right|},}$$

(10) $$avgPr\left(t\right)={\displaystyle \frac{1}{m\left(t\right)}\cdot \sum _{v=1}^{m\left(t\right)}p{r}_v\left(t\right),}$$

(11) $$avgRc\left(t\right)={\displaystyle \frac{1}{n}\cdot \sum _{v=1}^{n}r{c}_v\left(t\right).}$$

The optimal balance between precision and recall is calculated using the maximum F1-score, which for this study is calculated by calculating F1-scores with 0.01 step-sizes after transforming the Total Scores to a 0–1 range:

(12) $$F_{max}=\underset{t}{max}\left\{{\displaystyle \frac{2\cdot avgPr\left(t\right)\cdot avgRc\left(t\right)}{avgPr\left(t\right)+avgRc\left(t\right)}}\right\}.$$

The S_min score is calculated based on the remaining uncertainty ru and missing information mi, which is the IC of all false negatives and false positives, respectively. A threshold t step-size of 0.01 is used after transforming the Total Scores to a 0–1 range.

(13) $$S_{min}=\underset{t}{min}\sqrt{ru{\left(t\right)}^2+mi{\left(t\right)}^2,}$$

with ru as:

(14) $$ru\left(t\right)={\displaystyle \frac{1}{n}\sum _{v=1}^n{\sum }_{GO\in T_v-P_v\left(t\right)}IC\left(i\right),}$$

and mi as:

(15) $$mi\left(t\right)={\displaystyle \frac{1}{n}\sum _{v=1}^n{\sum }_{GO\in P_v\left(t\right)-T_v}IC\left(GO\right).}$$

Algorithm comparisons between Wei2GO and Argot

Because Wei2GO’s concept is based on that of Argot, it is important to highlight major algorithmic differences between the different software (Table 1). Both Argot2 and Argot2.5 supplement pfam2go annotations with the GOA GO terms of all proteins belonging to each Pfam entry, which Wei2GO does not. Wei2GO weighs GO annotations differently based on the associated evidence codes in the GOA database, which both Argot’s do not. Argot2.5 is the only software that performs taxonomy-based filtering of GO terms. And finally, Argot2 and Argot2.5 can only handle 5,000 and 10,000 proteins at a time, respectively.

To get an indication of the effect of the algorithmic differences between the software, Table 2 shows an overview of the correctly and incorrectly predicted protein-GO term pairs, irrespective of their prediction scores. The main difference between the software is the amount of false positive (FP) predictions. Argot2, which supplements pfam2go annotations, and to a lesser extend Argot2.5, which supplements pfam2go annotations and applies a taxonomic filter, show a substantially larger amount of FP predictions compared to Wei2GO whilst not showing a large increase in true positive (TP) predictions. Exceptions are Biological Processes (BPO) where Argot2.5 shows both a lower TP and FP count compared to Wei2GO, and Cellular Components (CCO) where Argot2 and Argot2.5 show an increase in TP count. This general trend suggests that supplementing Pfam annotations by adding GOA GO annotations from proteins belonging to the Pfam entry, which is not present in Wei2GO, adds quantity but not quality to the predictions.

Performance comparisons between Wei2GO and Argot

Several performance metrics are used to compare Wei2GO and Argot: precision-recall curves, F_max scores, S_min scores, and computational time. DeepGOPlus is added as a part of all these analyses. This deep learning-based and open source protein function prediction software is an excellent example of both predictive performance and usability. Adding this software in the assessments provides perspective as to any performance differences between Wei2GO and Argot. For the purpose of this analysis, all databases used by Wei2GO are the same versions used by Argot2 and Argot2.5. Each species was given to Argot2.5 as a separate input with a species level taxonomic code.

Precision recall curves were generated for biological process (BPO), molecular function (MFO), and cellular component (CCO) ontologies (Fig. 1). This analysis shows that Wei2GO is the best performer in both the BPO and MFO assessment, with most notably BPO having a substantial increase in precision. For CCO, Wei2GO is slightly outperformed by Argot2, while DeepGOPlus performs substantially better than both.

F_max scores were calculated based on the highest F1-score at any given threshold, and S_min scores were calculated based on the lowest S-score at any given threshold (Table 3). The F_max score is highest in Wei2GO for BPO, similar between Wei2GO and Argot2.5 for MFO, and is outperformed by both Argot2 and DeepGOPlus for CCO. Wei2GO shows the lowest S_min score for both BPO and MFO, and is outperformed by both Argot2 and DeepGOPlus for CCO.

Wei2GO was compared to Argot2.5 on the CAFA3 data using the F_max score, relative to the top three performers of the challenge (Table 4). Both Wei2GO and Argot2.5 are competitive, and show similar results to those in Table 3. One notable difference is a discrepancy in performance for the BPO category. Protein function annotation is an ‘open world’ problem, i.e. a GO term annotation can be incorrectly assigned as a false positive due to an as-of-yet unidentified function of a protein (Dessimoz, Škunca & Thomas, 2013). The test set used for Table 3 is more complete for BPO terms than the CAFA3 challenge set, with 6.5 average leaf terms per protein compared to 3.6 respectively, resulting in a larger F_max difference between Wei2GO and Argot2.5. These differences in GO terms per protein were marginal for MFO and CCO terms.

Computational times were compared in Table 5 based on the 4,067 proteins of the test set. The Argot2 and Argot2.5 web-servers were run without queue times, and all analyses were run on one CPU.

Implementation and availability

Wei2GO requires Python 3.6 or higher, and the NumPy package. Supplemental Files used to create the results in this paper, instructions on how to install the software, and instructions on how to run the software can be found on https://gitlab.com/mreijnders/Wei2GO. Additionally, Snakemake scripts are provided to run the full annotation pipeline, and update the databases to their latest versions.