# Directed acyclic coloured multigraphs for expert eukaryotic gene annotation

Sarah Djebali, Franck Delaplace, Hugues Roest Crollius

► To cite this version:

Sarah Djebali, Franck Delaplace, Hugues Roest Crollius. Directed acyclic coloured multigraphs for expert eukaryotic gene annotation. Journée Ouverte Biologie Informatique Mathématiques (JOBIM 2005), 2005, Paris, France. hal-00339795

HAL Id: hal-00339795

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

Submitted on 26 Aug 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished 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_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Directed acyclic coloured multigraphs for expert eukaryotic
gene annotation

Sarah Djebali¹,² Franck Delaplace² and Hugues Roest Crollius¹

¹ Laboratoire Dynamique et Organisation des Génomes, UMR CNRS 8541, École Normale Supérieure, 46 rue d’Ulm,
75005 PARIS, France
djebali@ens.fr, hrc@ens.fr

² Laboratoire de Méthodes Informatiques, UMR CNRS 8042, Université d'Evry Val d'Essonne, Genopole Tour Evry 2,
523 place des terrasses de l'Agora, F-91000 EVRY, France
delapla@lami.univ-evry.fr

**Abstract:** Automatic gene annotation is a basic requirement to identify the precise location of all the genes in a genome. Today however, the highest quality gene annotations remain those generated by human experts based on different reliable resources. Here we describe Exogean, a new automatic gene annotation method that combines two such resources, i.e. cross-species sequence alignments and mRNA sequences, using the same heuristic rules as human experts. Exogean is based on directed acyclic coloured multigraphs (DACMs) where nodes represent biological objects (exons, transcripts), and multiple edges between nodes represent relations between objects derived from human expertise. DACMs are read along paths that replicate the rules applied by human annotators when processing data. By representing relations and rules separately, the use of DACMs establishes a strict independence between the data and the heuristics applied to the data. Exogean is thus modular and flexible, and represents a generic framework that may capture any human expertise in the field of eukaryotic gene annotation. We describe the evaluation of Exogean on human chromosome 22 in comparison to eight other gene annotation methods, and we evaluate the influence of the evolutionary distance and the exhaustiveness of the input resources on these performances. The Exogean source code in OCaml is available on request at djebali@ens.fr.

# 1 Introduction

Ideally, the process of annotating protein coding genes (hereby referred to as “genes”) in a region of genomic DNA means locating the boundaries of all the genes it includes and for each, finding all the possible transcript variants. In practice, this procedure is very difficult in eukaryotic genomes for many reasons. First, eukaryotic genes are generally composed of a succession of exons and introns which makes their structure complex and highly variable. Second, genes cover only a small fraction of eukaryotic genomes (30% in mammals) and exons cover an even lower fraction (1 to 2% in mammals). Third, some eukaryotic genomes contain many pseudo-genes, i.e. non functional copies of genes often nested within genes and with similar compositions. Finally, each gene may give rise to many different transcripts depending on the spatial or temporal requirement for the corresponding protein function. Despite these difficulties, precise gene annotation is crucial for biomedical research: it is a basic requirement to link genotype and phenotypes in human and in model species and generally to focus the work of biologists and bio-informaticians on an essential functional part of the genome. Twenty eukaryotic genomes are now at least partially sequenced and assembled and each is commonly hundred millions or even billions nucleotides long: annotating genes in this massive amount of data requires a mathematical model.

Mathematical models have been proposed to automatically find genes in genomic DNA, either by similarity to expressed sequences, or by capturing our current biological understanding of genes in statistical algorithms.
The first are called *similarity* models and mostly consist in aligning transcribed sequences (ESTs, mRNAs) to a
genomic sequence in order to locate its genes (Est2genome [6], SIM4 [7], Spidey [17]), or to find all the possible
transcript variants of its genes in more sophisticated cases (ESTGenes [27], Multiassembly Problem [28]). The
second are called *ab initio* models and are usually divided in two categories [8,18]: the exon-based and the
signal-based approaches. The exon-based approach first finds and scores potential exons in a sequence and
then assembles them into a gene structure (GeneId [1], Fgene [3]), whereas the signal-based approaches, based
on Hidden Markov Models (HMM), find gene structures by considering them as successions of signals (i.e.
transcription and translation signals) separated by 'homogeneous' regions (Genscan [4], Genie [11], HMMgene
[5] and recently SNAP [32]). However pure similarity and *ab initio* methods are not sufficient to reliably
annotate a genomic sequence, mainly because the former critically depend on the quality and completeness
of the sequences used for comparison, while the latter remain too imprecise due to the lack of discriminative
power of signals, and thus predict too many incorrect genes.

Recently, with the emergence of comparative genomics, efforts have been made to combine the models
underlying the two approaches into a single framework. Among these are two categories: those that start from
similarity program outputs and build gene models on this basis using *ab initio* methods, and those that follow
an *ab initio* approach and try to incorporate results from similarity programs in the underlying model. In the
first category, the underlying *similarity* model can be based on protein sequences, such as in Genewise [33]
or on genomic sequences such as in SLAM [22], Doublescan [21], or Progen [12]. In the second category
the underlying *ab initio* model may be exon-based as in SGP [13] or signal-based as in Twinscan [14] or
Genomescan [15]. An altogether different approach is to consider the outputs of different methods (such as
those cited above) as resources and to combine them with the aim of obtaining better performances than either
method alone. In this more recent category are methods that only combine *ab initio* predictions as in Genescope
[9] and DIGIT [24], and methods that also incorporate similarity results such as in Gaze [19], Combiner [29],
Eugene'Hom [25], AIR [35] and Ensembl [30].

Despite these tremendous advances, automatic gene annotations are still considered as predictions that re-
quire validation by human experts, for instance when experimental work should be based upon them. In fact,
high quality reference gene annotations (e.g. those collated in the Vega repository [37]) are still considered
to be those generated by humans based on a number of resources (cross-species sequence alignments, mRNA
sequences, *ab initio* predictions, etc.). These resources are reviewed on a gene-by-gene basis using strict rules
based on a deep knowledge of both the data at hand and the biology associated with gene expression (trans-
scription, splicing, translation, etc). Based on these observations, we were interested in designing an automatic
annotation method that explicitly establishes the same relations between biological objects, and applies the
same rules upon them, as human experts do. The method must be able to deal with different resources and
different rules, both of which are susceptible to evolve with time. This is difficult since very often two heuristic
rules will involve two different data representations and processes. It is therefore instrumental to maintain a
strict independence between the determination of relationships between biological objects and the application
of heuristic rules to this data, although both define the expertise that we wish to derive from human knowl-
edge. To apply these principles, we have used Directed Acyclic Coloured Multigraphs in a new method called
Exogean (EXpert On GEne ANnotation).

DACMs applied to eukaryotic gene annotation

Although every human expert will approach the problem of annotating genes from slightly different angles
using their own specific knowledge and experience, we can define four general steps that humans follow when
confronted with a region of genomic DNA that must be annotated: (i) preprocessing of initial data to exclude
obvious artefacts, (ii) recurrent clustering and integration of data into gene models of increasing complexity,
(iii) identification of coding sequence (CDS) boundaries, (iv) final assessment and quality checks. It should be
noted that humans do not necessarily proceed along the four steps strictly in this order, and often jump back and forth between steps. The second step is the core of the annotation procedure, where heterogeneous and potentially conflicting data are integrated to generate one or several transcript models. This is performed in Exogeans by Directed Acyclic Coloured Multigraphs (DACMs), which integrate both the biological objects and the relations between them into a mathematical formalism. An overview of the Exogeans method is presented in figure 1:

<figure><img src="image_2.png" /><figcaption>**Figure 1.** Exogeans overview.</figcaption></figure>

- Preprocessing: normalises mRNA and protein alignments (hereafter called *HSPs*, following the BLAST terminology) by eliminating clear artefacts using heuristic rules (see supplementary data [39]).

- DACM expert annotation: for each gene, preprocessed mRNA and protein HSPs are integrated within one or several transcript models, in three steps. In each step, a graph is built where nodes are biological objects (exons or transcript) and edges between nodes are relations derived from human expertise. Each graph is then read along a path also defined by human expertise, resulting into a new graph where nodes are objects of increased complexity (see Methods). For each gene multiple alternative transcripts are found.

- CDS identification: a coding sequence (CDS) is defined for each transcript model using heuristic rules (see supplementary data [39]).

- Quality filter: transcript models with their associated CDS are classified into either final gene model corresponding to functional protein coding genes or into pseudogene models using heuristic rules (see supplementary data [39]).

While Exogeans is composed of four modules, this section will only cover the DACM algorithm since the other three are more intuitive serial applications of heuristic rules derived from basic biological expertise (see Supplementary Information [39]). In a DACM, vertices represent biological objects (protein or mRNA alignments,
proteins, mRNAs) and edges represent some known relations between these objects. Several edges may connect two vertices. A DACM is defined as follows:

**Definition (DACM)**

A Directed Acyclic Coloured Multigraph (DACM) is a triplet $\langle V, E, C \rangle$ where $V$ is a set of vertices, $C$ is a set of colours, $E \subseteq V \times V \times C$ is a set of directed edges, and such that there is no directed cycle in the graph.

Human experts do not directly build gene models from initial HSPs but rather tend to recurrently clusterise the data to form more and more complex objects. If we use DACMs to represent these objects and the relations between them, then in terms of graphs the expert protocol is equivalent to reducing the graphs. The following definition more formally describes a reduced graph.

**Definition (Reduced graph)**

Let $G = \langle V, E, C \rangle$ be a DACM, $G' = \langle V', E', C' \rangle$ is said to be the reduced graph of $G$ according to $\varphi: 2^V \to V'$ if and only if, let $P_G = \{V_1, \dots, V_n\}$ be a partition of $V^3$, then the two following properties are fulfilled:

- Each vertex of $V'$ is the result of the mapping by $\varphi$ of a set of vertices $V_i$:

$$ \forall v'_i \in V', \exists v_i \in P_G, \varphi(v_i) = v'_i $$

- Each edge of $E'$ is bound to the existence of an edge in $E$ with vertices that belong to two different sets of vertices $V_i$ and $V_j$:

$$ \forall (v'_i, v'_j, c) \in E', \exists v_i, v_j \in P_G, \exists (v_i, v_j) \in V_i \times V_j, (v_i, v_j, c) \in E. $$

Note that knowing a partition $P_G$ is equivalent to knowing a $\varphi$ function. It may happen that $P_G$ is not unique (for instance in the second DACM reduction (subsection 1.2)). In such cases all the reductions corresponding to the different partitions are applied, and each provides a different biological solution.

A subset of colours may also define a reduced graph, since any $C_1 \subseteq C$ uniquely defines $P_G = \{V_1, \dots, V_n\}$, where each $V_i$ is a maximal path in $G$ with all edges labeled with the colours of $C_1$. In the two first DACM reductions (subsections 1.1 and 1.2), we use this last definition of reduced graph and we associate it to the expression reduced graph according to a subset of colours.

As indicated on figure 1, the DACM expert annotation is composed of three main parts, each building more elaborate transcript models:

- the Single Molecule Clustering (or Level 1 Transcript Modelling),

- the Single Type Multi Molecule Clustering (or Level 2 Transcript Modelling),

- the Multi Type Multi Molecule Clustering (or Level 3 Transcript Modelling).

These three parts of the DACM expert annotation all use a DACM data structure and algorithm. They are described in more details in the next three subsections all divided into four parts: Introduction, DACM data structure, DACM algorithm and Example. In all the examples illustrating the DACM expert annotation circle nodes represent mRNAs and square nodes represent proteins. Note that the biological objects manipulated in the DACMs are ordered by beginning position on the genomic sequence.

## 1.1 Single Molecule Clustering

The Single Molecule Clustering (or Level 1 Transcript Modelling) assembles preprocessed HSPs from the same molecule into preliminary transcript models called level 1 transcript models (L1TM). HSPs from the same molecule and aligned to the same genomic sequence $S$, are not necessarily involved in the same gene on $S$, due to gene duplication and alignment artefacts. The single molecule clustering discards such situations by grouping preprocessed HSPs $h_1$ and $h_2$ from the same molecule only if they fulfil the following conditions:

3 A partition $\{V_1, \dots, V_n\}$ of $V$ fulfills the following properties:

$$ \cup_i V_i = V, i \neq j \Rightarrow V_i \cap V_j = \emptyset $$
- the genomic distance between $h_1$ and $h_2$ is less than a threshold (representing a maximal intron size),

- $h_1$ comes before $h_2$ on their common molecule sequence.

**DACM data structure** The single molecule clustering data structure is a DACM $\langle V, E, C \rangle$ where:

- $V$ is the set of preprocessed mRNA and protein HSPs,

- $C$ is the set of colours $\{\text{same-molecule}, \text{maximal-intron-size}, \text{genomic-molecule-order}\}$.

**DACM algorithm** The single molecule clustering algorithm reduces the graph according to $C$, i.e. to the *same-molecule*, the *maximal-intron-size* and the *genomic-molecule-order* colours.

**Example** Figure 2 illustrates the Single Molecule Clustering. Nodes represent preprocessed HSPs. Nodes $h_1$ to $h_4$, $h_5$ to $h_6$, $h_7$ to $h_9$, $h_{14}$ to $h_{16}$ and $h_{17}$ to $h_{18}$ are clustered in L1TM $r1$, $r2$, $r3$, $p3$ and $p_4$ respectively, since they represent maximal paths labeled with the three possible colours (*same-molecule*, *maximal-intron-size* and *genomic-molecule-order*). However the lack of an edge labeled with the *genomic-molecule-order* colour between $h_{11}$ and $h_{12}$ results in two L1TMs: $p_1$ (including $h_{10}$ and $h_{11}$), and $p_2$ (including $h_{12}$ and $h_{13}$).

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Single Molecule Clustering (or Level 1 Transcript Modelling).</figcaption></figure>

## 1.2 Single Type Multi Molecule Clustering

The *Single Type Multi Molecule Clustering* (or Level 2 Transcript Modelling), groups level 1 transcript models of the same molecule type (i.e. mRNA or protein), into more elaborate transcript models called level 2 transcript models (L2TM). Again mRNAs and proteins are treated differently by the algorithm: overlapping mRNA level 1 transcript models with a compatible splicing structure are grouped into one of potentially several alternative transcript classes, whereas for overlapping protein level 1 transcript models are grouped into a single transcript.

**DACM data structure** The single type multi molecule clustering underlying data structure is a DACM $G = \langle V, E, C \rangle$ where:

- $V$ is a set of mRNA and protein level 1 transcript models,

- $C$ is the set of colours $\{\text{genomic-overlap}, \text{extension}, \text{inclusion}\}$.
**DACM algorithm** For mRNAs the algorithm is inspired from the ClusterMerge algorithm used in ESTGenes for ESTs [27]: an alternative transcript class is built from all the mRNA level 1 transcript models included in a path from a root to a leaf following extension edges and including the “inclusion subgraph” of each vertex of the path. For proteins the algorithm is the reduction of the protein subgraph according to the *genomic-overlap* colour.

**Example** Figure 3 illustrates the Single Type Multi Molecule Clustering. Nodes represent L1TMs. mRNA nodes $r_1$, $r_2$ and $r_3$ are clustered in two L2TMs since they form two maximal mRNA paths labeled with the *extension* colour: $r_1$ to $r_3$, and $r_2$ to $r_3$. Protein nodes $p_1$ to $p_4$ are clustered in three L2TMs since they form three maximal protein paths labeled with the *genomic-overlap* colour: $p_1$ to $p_3$, $p_2$ and $p_4$.

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Single Type Multi Molecule Clustering (or Level 2 Transcript Modelling).</figcaption></figure>

## 1.3 Multi Type Multi Molecule Clustering

The *Multi Type Multi Molecule Clustering* (or Level 3 Transcript Modelling), groups level 2 transcript models into multimolecule transcript models. It is an important step since it merges the mRNA and the protein resources. More precisely it links every mRNA level 2 transcript model to all the protein level 2 transcript models of compatible structure that overlap with it on the genome. As stated before, mRNAs and proteins each have their own particular advantages which make their combination useful. In particular proteins provide the final transcript ORF, mRNAs provide the final transcript splicing structure and both resources may indicate a genomic extension (either in the 5' or the 3' direction).

**DACM data structure** The multi type multi molecule clustering underlying data structure is a bipartite DACM $G = \langle V, E, C \rangle$ where:

- $V$ is a set of mRNA and protein level 2 transcript models,

- $C$ is the set of colours {genomic-overlap, compatible-splicing-structure}.

**DACM algorithm** The multi type multi molecule clustering algorithm associates each mRNA level 2 transcript model of $G$ to the set of all its protein level 2 transcript model successors via the two colours of $C$, i.e. the *genomic-overlap* and the *compatible-splicing-structure* colours.

**Example** Figure 4 illustrates the Multi Type Multi Molecule Clustering. Nodes represent L2TMs. mRNA node $R_1$ is clustered with protein nodes $P_2$ and $P_3$ but not with protein node $P_1$, for $P_2$ and $P_3$ are successors of $R_1$ via two edges: one labeled with the *genomic-overlap* colour and one labeled with the *compatible-splicing-structure* colour, but $P_1$ is a successor of $R_1$ via only one edge labeled with the *genomic-overlap* colour. Following the same rules mRNA node $R_2$ is clustered with protein nodes $P_2$ and $P_3$.
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Multi Type Multi Molecule Clustering (or Level 3 Transcript Modelling).</figcaption></figure>

**1.4 Running Exogean**

Exogean is implemented in OCaml and has a time complexity in $\Theta(E)$, with E being the number of edges. Here, Exogean was tested on human chromosome 22 [10] from NCBI build 34. We aligned 6,791 Genbank [26] mRNA from human chromosome 22, keeping the alignment with the least mismatches in case of multiple hits, using BLAT [20] with the following parameters: -fine, -minIdentity=98, -trimHardA, -trimT, -ooc=0oc.11, -extendThroughN. We next aligned the entire set of 40,981 mouse proteins from the International Protein Index (IPI[38]; july 2004) using BLAT and the following parameters: -tileSize=7, -oneOff. Repeats from the genomic sequence of human chromosome 22 were masked using RepeatMasker. The combined computation time for the mRNA and protein alignments is 31 minutes on a single 2.4 GHz Xeon processor with 4 Gb memory. Exogean then computes annotations on chromosome 22 in approximately 3 minutes.

Chromosome 22 predicted annotations from Ensembl [30], Fgenesh++ [16], Acembly [42], SGP [13], Twinscan [14], Genscan [4] and Exoniphy [31] were downloaded from the UCSC website [40] (NCBI build 34). CCDS annotations from NCBI build 35 [41] were remapped on build 34 using the liftover utility from UCSC. Protein sets from other vertebrate species were downloaded from the International Protein Index (IPI) web site (february 2005).

**2 Results**

**2.1 Reference genes**

Evaluating the performances of Exogean requires a set of high quality reference gene models that have been reviewed by human experts. We chose protein coding genes from human chromosome 22 mostly because it is large enough to provide a diversity of human genes but is located on a small enough genomic segment so that each annotation (false or correct) can be examined individually. Human chromosome 22 (NCBI Build 34) is 34.8 Mb long and contains 406 genes (4176 exons, and excluding immunoglobulin genes) that possess an annotated coding sequence (CDS) free of stop codons, larger than 100 amino acids and multiple of three [23]. The sequence also contains 171 partial genes that did not meet the above criteria and 363 pseudogenes and immunoglobulin genes. These latter structures were not used as references, and any model overlapping either pseudogenes, partial genes or immunoglobulin genes was excluded from the analysis.

**2.2 Exogean Performances**

We have evaluated the performances of Exogean in comparison with a number of automatic gene annotation methods (Ensembl, SGP, Twinscan, Genscan, Acembly, Fgenesh++), one exon predicting tool (Exoniphy) and
one consensus set of high quality, semi-manually curated gene models (CCDS). To calculate performances in terms of specificity and sensitivity [2], a problem arises due to the fact that reference genes provide one transcript per gene while some methods (Exoge an, Ensembl, Acembly) may predict several alternative transcripts per gene. In addition, we wished to evaluate performances both in terms of genes predicted exactly (CDS start, all splice sites and stop) and in terms of predictions that simply overlap reference genes (boundaries overlap by at least one base). Therefore, we decided to consider all transcripts predicted for a given gene as counting for a single true positive when at least one transcript exactly matches (resp. overlap) the reference gene, and once as a false positive when all transcripts of a predicted gene do not match (resp. overlap) a reference gene. Sensitivity and specificity are thus expressed as:

$$Sn_{Exact} = \frac{Exact}{Exact+Overlap+False Negatives}$$

$$Sp_{Exact} = \frac{Exact}{Exact+Overlap+False Positives}$$

$$Sn_{Overlap} = \frac{Overlap}{Exact+Overlap+False Negatives}$$

$$Sp_{Overlap} = \frac{Overlap}{Exact+Overlap+False Positives}$$

Performances for all nine methods were calculated in the same way both for genes and for exons, except for Exoniphy where only exons were evaluated. Results are summarised in figure 5. The sensitivity of Exoge an

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Specificity and sensitivity of nine sources of gene annotations, in comparison with 406 protein coding genes on human chromosome 22. Methods are arranged from left to right by order of decreasing performances in terms of genes predicted exactly.</figcaption></figure>

is higher than for any other method, while the specificity is lower than for CCDS and slightly lower than for Ensembl predictions. The high specificity of CCDs is expected since it represents a set of consensus CDSs from several high quality sources including human expert curation (see methods).

In principle, a lower specificity implies a higher level of false positive predictions. On the other hand, the set of reference gene is by no means an exhaustive list of human chromosome 22 genes, especially as human mRNAs and mouse protein catalogues increase continually in quality and in quantity. We were thus interested to examine both Ensembl and Exoge an false positive annotations on a gene by gene basis to assess whether human curation could further classify them as new genes (i.e. absent from the set of reference genes) or true false predictions. Among the 47 false positive genes predicted by Exoge an, 30 were classified as true positives by human curation while 16 of the 39 false positives predicted by Ensembl were classified as true positives. Performances for Exoge an and Ensembl with this new data are presented in table 1.
<table><caption>**Table 1.** Exoge and Ensembl gene performances after manual curation.</caption><thead><tr><th rowspan="2"></th><th colspan="2">Exoge</th><th colspan="2">Ensembl</th></tr><tr><th>%</th><th>Nr of genes</th><th>%</th><th>Nr of genes</th></tr></thead><tbody><tr><td rowspan="2">Gene Overlap</td><td>Sp 95.94</td><td>-</td><td>94.28</td><td>-</td></tr><tr><td>Sn 92.20</td><td>402</td><td>89.81</td><td>379</td></tr><tr><td rowspan="2">Gene Correct</td><td>Sp 73.75</td><td>-</td><td>72.89</td><td>-</td></tr><tr><td>Sn 70.87</td><td>309</td><td>69.43</td><td>293</td></tr></tbody></table>

## 2.3 Impact of the resources

**Incremental data** The performances of Exoge are clearly dependent on the resources at hand in terms of mRNAs and proteins from different species. To examine the impact of the resources on the specificity and the sensitivity, we ran Exoge in multiple configurations. First, genes were predicted solely on the basis of incremental amounts of either randomly selected protein or randomly selected mRNA sequences (figure 6). Notably, the specificity of the annotations is not affected by the amount of data, indicating that Exoge

<figure><img src="image_7.png" /><figcaption>**Figure 6.** Variation in sensitivity and specificity of Exoge gene annotations that overlap reference genes. Incremental, randomly selected, non redundant fractions of the complete set of mRNA sequences (a) or protein sequences (b) were used.</figcaption></figure>

satisfactorily limits the noise inherent to large biological dataset. At the same time, the control of background noise is not performed at the cost of the detection of real genes, since sensitivity progresses as expected with increasing amounts of data. Interestingly, while human mRNA seem to reach saturation at approximately 70-80% of the current set of human mRNA mapped to chromosome 22, mouse proteins provide a linear increase in sensitivity up to the full set of sequences currently available. This suggest that further refinements and enrichments of the set of mouse predicted proteins is likely to continue to be informative for the annotation of human genes.

**Species Combination** When using protein sequences from a different species than human to predict genes in human genomic DNA, the evolutionary distance is expected to have a strong influence on the sensitivity of an annotation method such as Exoge. We evaluated the impact of the evolutionary distance using mouse, rat, chicken, Tetraodon, or zebrafish predicted protein sequences, as well as different combinations of these
protein sets. Results are shown in figure 7. Again specificity remains constant across the vertebrate spectrum

<figure><img src="image_8.png" /><figcaption>**Figure 7.** Performances with different combinations of species.</figcaption></figure>

tested here, which is expected since protein alignment conditions were calibrated on the closest representatives (rodents). Hence few spurious alignments are expected as evolutionary distance increases. Sensitivity decreases with increasing evolutionary distance from human, but clearly different sets of genes are missed by different species. For instance, a combination of chicken and Tetraodon proteins identifies more human genes (68.1 %) than chicken alone (62.5 %). However, the sensitivity provided by combining all five vertebrate predicted proteomes (82.7%) is not significantly better than that provided by mouse alone (82.5 %).

# 3 Discussion

We present here Exogean, a new eukaryotic gene annotation method that combines several resources (currently mRNA and protein alignments) using a multigraph based algorithm. The decisions taken by the algorithm during the merging of the data are explicitly based on human expertise as opposed to a statistical model. The use of DACMs allows the independence between the relations that link different biological objects, and the rules that are applied to these objects. Consequently, relations and rules can be changed, and new ones added, without modifying the basic DACM structure. More precisely, new edges (relationships between biological objects) can be added, upon which new rules for reducing graphs (merging data) can be applied. At a more practical level, many parameters such as thresholds for fusing HSPs or for defining CDSs are also easily configurable.

The performances of Exogean in its current configuration are similar to those of Ensembl, despite the fact that Exogean uses a more restricted source of data (human mRNAs and mouse proteins). However, Exogean is also fast (human chromosome 22 is annotated in about 30 mins including alignments; see methods) and portable. Exogean also generates all the details of the decision procedure thus providing, from the original evidences to the final gene model, a full traceability of the construction of the model.

Among the genes incorrectly identified by Exogean (false positives), most are due to the definition of the correct CDS within predicted transcripts. Exogean cannot currently predict partial genes, i.e. it will always start a CDS with a methionine, and thus truncate the beginning of the true CDS up to the first available methionine in cases when the 5' end of the CDS is not within the predicted transcript. Similarly, it will always try to find a CDS in a transcript, even if the latter only covers untranslated regions (UTRs). Often this CDS will be too short and will be rejected, but in cases where it is sufficiently long (e.g. longer than 100 amino acids, or
shorter but spanning an intron) it will be predicted as a CDS. Finally, a number of false positives are still due to pseudogenes being incorrectly annotated as protein coding genes. Although Exogean currently identifies 202 pseudogenes on chromosome 22 (of which 29% overlap reference pseudogenes), the rules used to define them need improving, for instance by taking into account the absence of transcript data, or the protein sequence conservation across stop codons and frameshifts.

Other methods that are specifically designed to combine protein and mRNA (and/or EST) sequence alignments to predict multiple transcripts per gene are Ensembl, AIR [35] and Eannot [34]. The Ensembl method is based on an annotation pipeline where protein alignments (obtained by Genewise [33]) are considered first, to which mRNA alignments (obtained by Exonerate [36]) are combined in order to potentially provide extensions. Eannot is a different pipeline approach where ESTs and mRNAs are clustered on the genome to define gene boundaries, and transcript models are then evaluated against information provided by proteins. AIR is closer to Exogean in the respect that a Directed Acyclic Graph (DAG) algorithm is used but where vertices are exons and edges are introns. ESTGenes [27], a fourth method that only considers EST sequences, also uses a DAG algorithm in a similar strategy to Exogean but again with a single edge between two vertices.

In contrast to pipelines, Exogean integrates both the data and the expertise on the data into the same graph representation, which provides greater flexibility when modifying the relations between objects and/or the rules applied to these objects. In addition, the data sources are considered on the same footing as opposed to a sequential (and inherently hierarchical) order, which means that more sources can easily be added if desired (e.g. ab initio predictions, genomic-genomic alignments). In comparison to conventional DAGs used in gene annotation the Exogean DACM represents, via the graph reduction mechanism, biological objects of increasing complexity in a way similar to the analysis performed by humans. This provides a greater intuitiveness when manipulating rules and relations that represent human expertise within the algorithm.

We therefore believe that Exogean represents a framework that can formally represents any human expertise in the field of gene annotation, and may thus be re-used by any laboratory with different experience and knowledge, or to predict genes in different species than human.

## Acknowledgements

We thank Jean-Marc Aury, Genoscope for assistance with computing the program performances.

## References

[1] R. Guigo, S. Knudsen, N. Drake, T. Smith, Prediction of gene structure. *J. Mol. Biol.*, 226:141-157, 1992.

[2] M. Burset, R. Guigo, Evaluation of gene structure prediction programs. *Genomics*, 34:353-367, 1996.

[3] V. Solovyev, A. Salamov, The Gene-Finder computer tools for analysis of human and model organisms genome sequences. *AAAI ISMB1997*, 294-302, 1997.

[4] C. Burge, S. Karlin, Prediction of complete gene structures in human genomic DNA. *J. Mol. Biol.*, 268:78-94, 1997.

[5] A. Krogh, Two methods for improving performance of an HMM and their application for gene finding. *AAAI ISMB1997*, 179-186, 1997.

[6] R. Mott, EST-GENOME: a program to align spliced DNA sequences to unspliced genomic DNA. *Comput. Appl. Biosci.*, 13:477-478, 1997.

[7] L. Florea, G. Hartzell, Z. Zhang, G.M. Rubin, W. Miller, A computer program for aligning a cDNA sequence with a genomic DNA sequence. *Gen. Res.*, 8:967-974, 1998.

[8] R. Guigo, Assembling genes from predicted exons in linear time with dynamic programming. *J. Comput. Biol.*, 5:681-702, 1998.

[9] K. Murakami, T. Takagi, Gene recognition by combination of several gene-finding programs. *Bioinformatics*, 14:665-675, 1998.

[10] I. Dunham, N.Shimizu, B.A. Roe, S. Chissoe, A.R. Hunt, J.E. Collins, R. Bruskiewich, D.M. Beare, M. Clamp, L.J. Smink, et al., The DNA sequence of human chromosome 22. *Nature*, 402:489-95, 1999.
[11] M.G. Reese, D. Kulp, H. Tammana, D. Haussler, Genie-Gene finding in Drosophila melanogaster. *Gen. Res.*, 10:529-538, 2000.

[12] P.S. Novichkov, M.S. Gelfand, A.A. Mironov, Gene recognition in eukaryotic DNA by comparison of genomic sequences. *Bioinformatics*, 17:1011-1018, 2001.

[13] T. Wiehe, S. Gebauer-Jung, T. Mitchell-Olds, R. Guigo, SGP-1: prediction and validation of homologous genes based on sequence alignments. *Gen. Res.*, 11:1574-1583, 2001.

[14] I. Korf, P. Flicek, D. Duan, M.R. Brent, Integrating genomic homology into gene structure prediction. *Bioinformatics*, 17:S140-S148, 2001.

[15] R-F Yeh, L.P. Lim, C.B. Burge, Computational inference of homologous gene structures in the human genome. *Gen. Res.*, 11:803-816, 2001.

[16] V.V. Solovyev, Statistical approaches in Eukaryotic gene prediction. *the Handbook of Statistical Genetics, John Wiley & Sons, Ltd.*, 83-127, 2001.

[17] S.J. Wheelan, D.M. Church, J.M. Ostell, Spidey: a tool for mRNA-to-Genomic alignments. *Gen. Res.*, 1952-1957, 2001.

[18] C. Mathe, M-F. Sagot, T. Schiex, P. Rouze, Current methods of gene prediction, their strengths and weaknesses. *Nucl. Acids Res.*, 30:4103-4117, 2002.

[19] K.L. Howe, T. Chothia, R. Durbin, GAZE: a generic framework for the integration of gene-prediction data by dynamic programming. *Gen. Res.*, 12:1418-1427, 2002.

[20] J.W. Kent, BLAT—The BLAST-Like Alignment Tool. *Gen. Res.*, 12:656-664, 2002.

[21] I.M. Meyer, R. Durbin, Comparative ab initio prediction of gene structures using pair HMMs. *Bioinformatics*, 18:1309-1318, 2002.

[22] M. Alexandersson, S. Cawley, L. Pachter, SLAM: Cross-species Gene finding and Alignment with a Generalized Pair Hidden Markov Model. *Gen. Res.*, 13:496-502, 2003.

[23] J.E. Collins, M.E. Goward, C.G. Cole, L.J. Smink, E.J. Huckle, S. Knowles, J.M. Bye, D.M. Beare, I. Dunham, *Gen. Res.*, 13:27-36, 2003.

[24] T. Yada, T. Takagi, DIGIT: a novel gene finding program by combining gene-finders. PSB2003, 2003.

[25] S. Foissac, P. Bardou, A. Moisan, M-J. Cros, T. Schiex, Eugene'Hom: a generic similarity-based gene finder using multiple homologous sequences. *Nucl. Acids Res.*, 31:3742-3745, 2003.

[26] D.A. Benson, I. Karsch-Mizrachi, D.J. Lipman, J. Ostell, D.L. Wheeler, GenBank. *Nucl. Acids Res.*, 31:23-27, 2003.

[27] E. Eyras, M. Caccamo, V. Curwen, M. Clamp, ESTGenes: Alternative splicing from ESTs in Ensembl. *Gen. Res.*, 14:976-987, 2004.

[28] Y. Xing, A. Resch, C. Lee, The multiassembly problem: reconstructing multiple transcript isoforms from EST fragment mixture. *Gen. Res.*, 14:1-16, 2004.

[29] J.E. Allen, M. Pertea, S.L. Salzberg, Computational gene prediction using multiple sources of evidence. *Gen. Res.*, 14:142-148, 2004.

[30] V. Curwen, E. Eyras, T.D. Andrews, L. Clarke, E. Mongin, S.M.J. Searle, M. Clamp, The Ensembl automatic gene annotation system. *Gen. Res.*, 14:942-950, 2004.

[31] A. Siepel, D. Haussler, Computational identification of evolutionarily conserved exons. RECOMB '04, 2004.

[32] I. Korf, Gene finding in novel genomes. BMC Bioinformatics, 5:59-56, 2004.

[33] E. Birney, M. Clamp, R. Durbin, Genewise and Genomewise. *Gen. Res.*, 14:988-995, 2004.

[34] L. Ding, A. Sabo, N. Berkowicz, R.R. Meyer, Y. Shotland, M.R. Johnson, K.H. Pepin, R.K. Wilson, J. Spieth, EAnnot: A genome annotation tool using experimental evidence. *Gen. Res.*, 14:2503-2509, 2004.

[35] L. Florea, V. Di Francesco, J. Miller, R. Turner, A. Yao, M. Harris, B. Walenz, C. Mobarry, G.V. Merkulov, R. Charlab, I. Dew, Z. Deng, S. Istrail, P. Li, G. Sutton, Gene and alternative splicing annotation with AIR. *Gen. Res.*, 15:54-66, 2005.

[36] G.S. Slater, E. Birney, Automated generation of heuristics for biological sequence comparison. BMC Bioinformatics, 6:31, 2005.

[37] J.L. Ashurst, C.-K. Chen, J.G.R. Gilbert, K. Jekosh, S. Keenan, P. Meidl, S.M. Searle, J. Stalker, R. Storey, S. Trevanion, The vertebrate genome annotation (Vega) database. *Nucl. Acids Res.*, 33:D459-D465, 2005.

[38] P.J. Kersey, J. Duarte, R. Apweiler, IPI - International Protein Index. *Nucl. Acids Res. Molecular Biology Database Collection*, 696, 2005.

[39] http://www.biologie.ens.fr/lgmldog/exogean/exogeansuppldata.html

[40] http://www.genome.ucsc.edu/.

[41] http://www.ncbi.nlm.nih.gov/CCDS/.

[42] http://www.ncbi.nlm.nih.gov/IEB/Research/Acembly.