e79 Nucleic Acids Research, 2009, Vol. 37, No. 11 hypergeometric distribution and the probability of seeing x red balls simply by chance is given by P x m1 x m1 x n m2 x m2 n x x x Pr½m1 ¼ xjm1 þ m2 ¼ n ¼ where maxð0, n m2Þ x minðm1, nÞ. The odds ratio, , is a measure of bias and is equivalent to the probability ratio of red over white balls: p1ð1 p2Þ ¼ p2ð1 p1Þ We reasoned that when the amount of annotation asso- ciated with a gene list is higher than expected by chance, the sampling will be biased in favour of certain types of tokens. As a consequence, some common words or non- speciﬁc terms shared by most abstracts, such as ‘cell’, ‘ana- lysis’ and ‘molecule’, will have higher probabilities of being selected along with genes in the gene list. As such, token signiﬁcance inferred from the classical hypergeo- metric distribution is likely to be misleading and errone- ous because this approach does not account for the excess annotation. Here, we propose a solution to this problem, which is based on the extended hyper-geometric distribu- tion model: (1) Odds ratio estimation: Let the number of genes anno- tated by token T in the background be m1, and the number of genes associated with tokens other than T in the background be m2. Given a gene list of size n, the odds ratio of T can be calculated empirically by: ~ yðm2 n þ yÞ ¼ ðm1 yÞðn yÞ where y is the mean number of genes we expect to see in the gene list just by chance. We determine y by ﬁtting a polynomial regression line through the token frequencies data: x1 þ x2 þ x3 þ x4 þ x5 þ x6 þ x7 y In this model, the explanatory variable x represents Chip frequency (equivalent to m1); the response vari- able y represents List frequency. The best-ﬁtting curve is obtained by the method of least squares implemented in the R function lm. For each token, the ﬁtted value of List frequency expected for a given Chip frequency is determined from the best-ﬁtting curve. This ﬁtted value is a good approximation for y and is used to calculate the odds ratio. (2) P-value calculation: Once the odds ratio has been determined, we use the pFNCHypergeo function implemented in the BiasedUrn R package to calcu- late a P-value for each token based on the extended hypergeometric distribution. Consensus gene age computation In order to gain an approximate measure of the ‘age’ (i.e. the length of time a gene has been known and PAGE 4 OF 11 researched upon) for genes on the HG-U133A array we used a combination of PubMed, OMIM (Online Mendelian Inheritance in Man) and HGNC (HUGO Gene Nomenclature Committee) to collate the following information for each gene: (1) the date of the earliest PubMed article that described the gene, (2) the date of the earliest article cited in an OMIM record for which the gene is described, (3) the date on which the gene ﬁrst appeared in an OMIM record and (4) the date on which the gene ﬁrst approved by HGNC. Based on these dates, we derived two ‘age estimates’ for each gene: literature- and database-based age. Literature-based age was calculated by averaging the values from (1) and (2), and it represents the date on which a gene was ﬁrst referred to in the scientiﬁc literature. Database-based age was calculated as the mean of (3) and (4), and represents the date when a gene symbol was ﬁrst approved or integrated into the public databases. There is a good correlation between the literature- and database- based ages (data not shown); therefore we calculated a ﬁnal consensus age for each gene by taking average of the literature- and database-based ages. Platform and availability The algorithms for Outlier and ExtendedHG have been implemented in a public web server PAKORA (http:// www.pakora.cf.ac.uk/pakora.php) for interactive use. The source codes (written in R and Perl) are available from the author upon request. A summary of the litera- ture-derived gene lists that we have used in this project can be found in Supplementary Data 2. These gene lists are available for download from PAKORA. The R scripts used in this analysis were developed and tested under R-2.6 and BioConductor-2.1. The Perl scripts were developed and tested under Perl v5.8.7. All analysis was performed on a Windows PC with a 2.8 GHz processor and 2 GB of RAM. RESULTS We began by analysing a list of interferon-stimulated genes (ISG) derived from the biological data reported in Sanda et al. (6), using the classical hypergeometric distri- bution-based ORA approach. The ISG gene list was used as a testbed because it constitutes a relatively simple and well-studied example of transcriptional regulation. PubMed articles associated with genes represented on the Aﬀymetrix HG-U133A array were collected and ﬁl- tered to give a text corpus consisting of 107 517 abstracts (70% of the unﬁltered collection). These abstracts were tokenised and stemmed to produce 220 290 unique single-word tokens for mining. Of these, 9486 tokens are associated with the ISG gene list. Classical hypergeometric test-based ORA is affected by annotation bias Initial use of the classical hypergeometric distribution- based method produced encouraging results when applied to the ISG gene list. In total, 94 tokens were called sig- niﬁcantly enriched after Bonferroni correction (Table 1).