CoXpress: differential co-expression in gene expression data
© Watson; licensee BioMed Central Ltd. 2006
Received: 21 July 2006
Accepted: 20 November 2006
Published: 20 November 2006
Traditional methods of analysing gene expression data often include a statistical test to find differentially expressed genes, or use of a clustering algorithm to find groups of genes that behave similarly across a dataset. However, these methods may miss groups of genes which form differential co-expression patterns under different subsets of experimental conditions. Here we describe coXpress, an R package that allows researchers to identify groups of genes that are differentially co-expressed.
We have developed coXpress as a means of identifying groups of genes that are differentially co-expressed. The utility of coXpress is demonstrated using two publicly available microarray datasets. Our software identifies several groups of genes that are highly correlated under one set of biologically related experiments, but which show little or no correlation in a second set of experiments. The software uses a re-sampling method to calculate a p-value for each group, and provides several methods for the visualisation of differentially co-expressed genes.
coXpress can be used to find groups of genes that display differential co-expression patterns in microarray datasets.
Microarrays have become a standard tool for the exploration of global gene expression changes at the cellular level . Data analysis often includes the use of a statistical test, such as a t-test or analysis of variance, to find genes differentially expressed in one set of conditions when compared to another, or the use of clustering algorithms in order to find groups of genes which behave similarly over a number of experiments . However, these techniques may not detect differential co-expression patterns that exist between two biological states.
Statistical tests, such as the t-test or ANOVA, identify genes that are differentially expressed under one or more conditions. The output of such tests is a simple list of genes, with an associated test statistic and p-value . There is no indication of which genes may be interacting with one another. Alternatively, clustering algorithms are often used to find groups of genes which display similar expression profiles across a dataset, and these clusters are subsequently analyzed visually for patterns of interest [4, 5]. Eisen et al used hierarchical cluster analysis to determine groups of co-expressed genes, and found that genes within those groups were functionally related , and the use of hierarchical cluster analysis is now a standard technique for analysing microarray data [2, 6]. Yeung et al  assessed the use of hierarchical clustering to find groups of co-regulated genes. Various clustering algorithms were used on a number of datasets, and the results evaluated by determining those genes that share a common transcription factor. Of the algorithms tested, MCLUST  and two hierarchical methods (based on the pearson correlation coefficient) showed the highest coincidence of correlated and co-regulated genes.
However, genes which show highly correlated patterns of expression in one biological state, but not in another, may not be highly correlated across the entire dataset, and therefore would not be associated with one another if a clustering algorithm is used. Variation may exist in the expression of a gene in different groups of individuals due to the presence of sub-populations, and this may lead to that gene being grouped incorrectly. Furthermore, clustering algorithms do not provide methods to identify groups that are behaving differently in different biological conditions.
Recent work has concentrated on alternative approaches to the discovery of co-expressed genes. Li  describes a method whereby genes whose expression is associated with differential co-expression patterns in other pairs of genes may be discovered, and Lai et al  describe a conceptually similar method whereby pairs of genes that display differential co-expression patterns between the normal and cancerous state may be discovered. Other approaches have centred on the construction of large gene co-expression networks. Lee et al  analysed 60 human microarray data sets to construct gene co-expression networks conserved across multiple data sets, and Stuart et al  constructed a gene co-expression network across different organisms, indicating that such relationships are evolutionarily conserved. However, neither of the above attempted to find group of genes differentially co-expressed between different conditions. Choi et al  tackled this problem by constructing normal and tumour co-expression networks from a variety of public datasets, comparing the results to find differences in co-expression patterns associated with cancer. In all of these cases, the networks were built by comparing genes pairwise, using some variation of the pearson correlation coefficient, to determine if a co-expression link exists between the two genes. These links were then joined to form a co-expression network.
Cluster analysis and network construction can be thought of as alternative methods for finding co-expressed genes. Whereas networks concentrate on conserved, pairwise comparisons, there is no guarantee that genes that are close in the network, but are not directly linked, have correlated expression profiles. Alternatively, cluster analysis produces groups of genes that are correlated above a certain level, defined by where the tree is cut and the clustering algorithm, but there is no indication of which particular pairs of genes are interacting. Kostka and Spang  described the first method to investigate differentially co-expressed groups of genes, using an additive model for scoring gene-gene co-expression and then a stochastic search algorithm to find groups of genes showing differential co-expression patterns. Jen et al  have produced ACT, the Arabidopsis Co-expression Tool, which allows users to calculate co-expression across user-defined data sets and uses a correlation cut-off to define groups of genes.
Here we describe coXpress, a simple and easy to use package that allows users to explore differential co-expression in an intuitive way. The package is aimed at biologists who want to analyse differential co-expression in their data set, which can be achieved in just 4 simple commands once the data has been loaded. CoXpress uses hierarchical cluster analysis to explore the relationship between genes, cutting the tree to form groups of genes that are co-expressed. This is an intuitive approach that many biologists are familiar with. CoXpress then uses a resampling approach to find those groups that are co-expressed in one set of experiments and not in another. The package should be used as first step in the analysis of co-expression, and is designed to complement the approaches described above.
CoXpress is released as a package for R. R is a freely available, open-source statistical package  that is widely used in the biological community. R has very powerful statistical and graphical capabilities, and many add-on packages are freely available. The bioconductor project [17, 18] provides a huge number of add-on packages for R, covering a wide range of biological data analysis applications, and the implementation of coXpress in R provides seamless integration with many of these packages. CoXpress is written in the native R language and has been fully tested on both windows and linux. R is available for windows, linux, unix and MacOS (including MacOS X).
The input for coXpress is a matrix of data, with rows representing genes and columns representing microarrays. The R data.frame object is most convenient, and can be created by reading in a text file (using the read.table function), an Excel spreadsheet (using the RODBC library) or from existing R objects, created by the packages from the bioconductor project such as affy , limma  or marray .
The genes are first clustered based on their expression values in a subset of experiments (termed subset 1), using the cluster.gene function. This function wraps the dist, cor and hclust functions that are built in to R, and thus provides a simple interface to hierarchical clustering. When a correlation coefficient is used as the distance measure, the distance measure is calculated as 1 - r, where r is the pearson correlation coefficient. The resulting tree is cut at a user-defined value, using the cutree function, to form groups of genes that are co-expressed in subset 1. These groups are then examined in both subset 1 and a second set of experiments, defined by the user, which we will term subset 2.
Groups of size 1 are ignored as there can be no co-expression. Groups of size two are handled by the cox.pairs function. The cox.pairs function uses the cor.test function in R to test if the genes are significantly correlated in subset 1 and subset 2. Thus, a pair of genes significantly correlated in subset 1 and not significantly correlated in subset 2 can be described as differentially co-expressed.
1 ≤ r ≤ +1
where 1 represents positive correlation, -1 represents negative correlation and 0 represents no relationship. The t-statistic is used here to summarise the "difference from zero" of a group of pairwise correlation coefficients. The exact formula for this is:
where x is the vector of unique, pairwise correlation coefficients, is the mean of x and se(x) denotes the standard error of x. A group of highly correlated genes will have a mean correlation close to 1 and a small standard error, resulting in a large value for t. However, a group of uncorrelated genes will have a mean close to 0 and a relatively large standard error, resulting in a small value for t. The observed t statistic is compared against m random t statistics in order to calculate a probability of randomness.
The AML/ALL leukaemia dataset
The utility of coXpress is demonstrated using gene expression data from the leukaemia microarray study of Golub et al . This dataset represents gene expression measurements from 38 tumour mRNA samples, 27 acute lymphoblastic leukaemia (ALL) cases and 11 acute myeloid leukaemia (AML) cases. The HU6800 Affymetrix array was used, which contains 6800 probesets. The dataset has been filtered such that genes with negative values in any sample have been removed, resulting in 2568 genes.
Using coXpress, the genes were first clustered according to their expression levels in the 27 ALL samples, using the cluster.gene function. The distance measure used was 1 - r, where r is the pearson correlation coefficient. The resulting tree was cut at a distance of 0.4, representing a correlation coefficient of 0.6, using the cutree function.
These groups were then examined in both the ALL and AML cases using the coXpress function. The observed t statistics in all cases were compared with the t statistics generated by randomly resampling the dataset 10,000 times for each group size. The resulting table contains one row for each group.
To test the robustness of the method to outliers, a bootstrapping approach was used. Each group was re-tested 1000 times, each time randomly selecting 75% of the observations for each leukaemia subtype (20 AML cases and 8 AML cases). The number of times each group was found to be differentially co-expressed by the coXpress method was recorded.
Differentially co-expressed groups from the Golub dataset
Most changed pairwise correlation coefficients between the ALL and AML subsets in group 62
Cor Group 1
Cor Group 2
Least changed pairwise correlation coefficients between the ALL and AML subsets in group 62
Cor Group 1
Cor Group 2
The GOHyperG function of the GOstats package  was used to find GO terms over-represented in the differentially co-expressed groups. Group 3, with 34 members, is enriched for GO terms for lymph node development, cell organisation and biogenesis, and protein biosynthesis and transport. Group 62, which has 7 members, is enriched for GO terms for methyltransferase activity, DNA modification, protein transport and DNA and protein methylation. Group 121, with 11 members, is enriched for GO terms for nucleotidase activity, and RNA splicing, processing and metabolism.
The ALL subtype dataset
This dataset is from the Acute Lymphoblastic Leukaemia study by Yeoh et al . Six subtypes of ALL leukaemias are represented in 248 cases. The six subtypes are T-ALL, E2A-PBX1, BCR-ABL, TEL-AML1, MLL rearrangement, and hyperdiploid >50. The HG_U95Av2 Affymetrix microarray was used which contains 12,600 probesets. The dataset has been filtered such that genes with negative values in any sample have been removed, resulting in 1516 genes present in the dataset.
Using coXpress, the genes were first clustered according to their expression levels in the BCR-ABL samples, using the cluster.gene function. The distance measure used was 1 - r, where r is the pearson correlation coefficient. The resulting tree was cut at a distance of 0.5, representing a correlation coefficient of 0.5, using the cutree function. These groups were then examined in both the BCR-ABL and T-ALL subsets.
Differentially co-expressed pairs in the ALL subtype dataset
Groups of N ≥ 3 were analysed in the BCR-ABL and T-ALL subsets using the coXpress function. The observed t statistics in all cases were compared with the t statistics generated by randomly resampling the dataset 10,000 times for each group size. The resulting table contains one row for each group.
To test the robustness of the method to outliers, a bootstrapping approach was used. Each group was re-tested 1000 times, each time randomly selecting 75% of the observations for each leukaemia subtype. The number of times each group was found to be differentially co-expressed by the coXpress method was recorded.
Differentially co-expressed groups from the ALL subtype dataset
The GOHyperG function of the GOstats package  was used to find GO terms over-represented in the differentially co-expressed groups. Group 47 with 16 members, is enriched for GO terms for hormone catabolism, glucocorticoid receptor signalling and glucocorticoid catabolism. Group 31 with 10 members contains two probes for a gene in the RAS oncogene family, and is enriched for GO terms for oxidoreductase activity and ubiquitin activating enzyme activity. Finally, group 89 with 13 members contains genes annotated as B-cell lymphoma and cancer susceptibility genes, as well as genes enriched for GO terms for endothelial cell migration, regulation of cell motility and migration, angiostatin binding and regulation of blood vessel endothelial cell migration.
It is clear that coXpress is capable of finding differentially co-expressed groups of genes in both data sets. The groups presented above are extremely highly correlated in one subset of experiments, yet show little correlation in another subset. Furthermore, these patterns of correlation are shown to be non-random in the first subset, and no different from random in the second subset. The results show that it is the overall correlation structures of these groups that have changed significantly and some pairs of genes are still highly correlated in the second subset. It is important that each group is examined using the inspect.groups function in order to determine which of the pairs of genes are still correlated and which are not. The groups found by coXpress could also feed into the network construction technique described by Choi et al  to determine which pairwise relationships are conserved and which are not. One would expect the differences between ALL and AML leukaemia in the Golub dataset to be larger than those between different ALL subtypes in the Yeoh dataset, and the fact that coXpress can still find groups with such different correlation structures demonstrates the power of the method.
The use of hierarchical cluster analysis, followed by cutting the tree, is an intuitive approach and one that is familiar to biologists. However, it has limitations. For example, each gene may only be in one group, which does not ring true for biological systems, where many genes have multiple functions. Also, the choice of where to cut the tree is arbitrary. A high cut-off will produce many small groups of genes that are very highly correlated, whereas a lower cut-off will produce fewer groups, of larger size, which are not as highly correlated. In reality the user must use a range of different cut-offs to see which performs best with their dataset. Other clustering algorithms, such as MCLUST [7, 8], have been shown to out-perform hierarchical cluster analysis, however, there is no reason why these algorithms could not be used to define the groups of genes prior to running the coXpress function.
There are several directions in which the software can be developed. At present, the user defines which subsets of experiments are analysed, however it is possible that coXpress could suggest, or improve, these groupings using an approach such as random forests or genetic algorithms. This may allow researchers to discover sub-populations in the system under study. Integration of other clustering algorithms with coXpress, such as MCLUST, may also improve the performance of the software. In particular, clustering or grouping algorithms that allow genes to be present in more than one group may be advantageous. Finally, the integration of network construction algorithms would allow researchers to further analyse and visualise the differentially co-expressed groups discovered by coXpress.
We describe coXpress, an open-source R package that allows researchers to analyse differential co-expression patterns in DNA microarray data. CoXpress contains several methods for the discovery and visualisation of differentially co-expressed genes. We show how coXpress can be used to find groups of differentially co-expressed genes in two publicly available microarray datasets. The groups found are shown to be highly correlated in one subset of experiments, yet show little or no correlation in a second subset of experiments. A comparison against random distributions is used to obtain a p-value for the co-expression of the genes in different subsets.
Availability and requirements
Project Name: coXpress
Project Home Page: http://coxpress.sf.net
Operating Systems: Windows, Linux
Programming Language: R
Other Requirements: R, gplots, gtools, gdata (for heatmaps), hu6800, hgu95av2, plotrix (for examples)
License: GNU GPL
acute lymphoblastic leukaemia
acute myeloid leukemia
T lineage leukaemias
B lineage leukemias that contain t(1;19)
B lineage leukemias that contain t(9;22)
B lineage leukemias that contain t(12;21)
- MLL rearrangement:
B lineage leukemias that contain rearrangements in the MLL gene on chromosome 11, band q23
This work was funded by the core strategic grant of the Institute for Animal Health, provided by the BBSRC.
- Schena M, Shalon D, Davis RW, Brown PO: Quantitative monitoring of gene expression patterns with a complementary DNA microarray. Science 1995, 270: 467–470. 10.1126/science.270.5235.467View ArticlePubMedGoogle Scholar
- Slonim D: From patterns to pathways: gene expression data analysis comes of age. Nat Genet 2002, 32(Suppl):502–508. 10.1038/ng1033View ArticlePubMedGoogle Scholar
- Stekel D: Microarray Bioinformatics. Cambridge University Press; 2003.View ArticleGoogle Scholar
- Eisen MB, Spellman PT, Brown PO, Botstein D: Cluster analysis and display of genome-wide expression patterns. Proc Natl Acad Sci USA 1998, 95: 14863–14868. 10.1073/pnas.95.25.14863PubMed CentralView ArticlePubMedGoogle Scholar
- Spellman PT, Sherlock G, Zhang MQ, Iyer VR, Anders K, Eisen MB, Brown PO, Botstein D, Futcher B: Comprehensive identification of cell cycle-regulated genes of the yeast Saccharomyces cerevisiae by microarray hybridization. Mol Biol Cell 1998, 9: 3273–3297.PubMed CentralView ArticlePubMedGoogle Scholar
- Quackenbush J: Computational analysis of microarray data. Nat Rev Genet 2001, 2(6):318–27. 10.1038/35076576View ArticleGoogle Scholar
- Yeung KY, Medvedovic M, Bumgarner RE: From co-expression to co-regulation: how many microarray experiments do we need? Genome Biol 2004, 5(7):R48. 10.1186/gb-2004-5-7-r48PubMed CentralView ArticlePubMedGoogle Scholar
- Yeung KY, Fraley C, Murua A, Raftery AE, Ruzzo WL: Model-based clustering and data transformations for gene expression data. Bioinformatics 2001, 17(10):977–87. 10.1093/bioinformatics/17.10.977View ArticlePubMedGoogle Scholar
- Li KC: Genome-wide coexpression dynamics: Theory and application. Proc Nat Acad Sci USA 2002, 99(26):16875–80. 10.1073/pnas.252466999PubMed CentralView ArticlePubMedGoogle Scholar
- Lai Y, Wu B, Chen L, Zhao H: A statistical method for identifying differential gene-gene co-expression patterns. Bioinformatics 2004, 20(17):3146–55. 10.1093/bioinformatics/bth379View ArticlePubMedGoogle Scholar
- Lee HK, Hsu AK, Sajdak J, Qin J, Pavlidis P: Coexpression analysis of human genes across many microarray data sets. Genome Res 2004, 14(6):1085–94. 10.1101/gr.1910904PubMed CentralView ArticlePubMedGoogle Scholar
- Stuart JM, Segal E, Koller D, Kim SK: A gene-coexpression network for global discovery of conserved genetic modules. Science 2003, 302(5643):249–55. 10.1126/science.1087447View ArticlePubMedGoogle Scholar
- Choi JK, Yu U, Yoo OJ, Kim S: Differential coexpression analysis using microarray data and its application to human cancer. Bioinformatics 2005, 21(24):4348–55. 10.1093/bioinformatics/bti722View ArticlePubMedGoogle Scholar
- Kostka D, Spang R: Finding disease specific alterations in the co-expression of genes. Bioinformatics 2004, 20(Suppl 1):I194-I199. 10.1093/bioinformatics/bth909View ArticlePubMedGoogle Scholar
- Jen CH, Manfield IW, Michalopoulos I, Pinney JW, Willats WGT, Gilmartin PM, Westhead DR: The Arabidopsis co-expression tool (ACT): a WWW-based tool and database for microarray-based gene expression analysis. The Plant Journal 2006, 46: 336–348. 10.1111/j.1365-313X.2006.02681.xView ArticlePubMedGoogle Scholar
- Gentleman RC, Carey VJ, Bates DM, Bolstad B, Dettling M, Dudoit S, Ellis B, Gautier L, Ge Y, Gentry J, Hornik K, Hothorn T, Huber W, Iacus S, Irizarry R, Leisch F, Li C, Maechler M, Rossini AJ, Sawitzki G, Smith C, Smyth G, Tierney L, Yang JY, Zhang J: Bioconductor: open software development for computational biology and bioinformatics. Genome Biol 2004, 5(10):R80. 10.1186/gb-2004-5-10-r80PubMed CentralView ArticlePubMedGoogle Scholar
- Gautier L, Cope L, Bolstad BM, Irizarry RA: affy: analysis of Affymetrix GeneChip data at the probe level. Bioinformatics 2004, 20(3):307–15. 10.1093/bioinformatics/btg405View ArticlePubMedGoogle Scholar
- Smyth GK: Linear models and empirical bayes methods for assessing differential expression in microarray experiments. Stat Appl Genet Mol Biol 2004, 3(1):Article3.PubMedGoogle Scholar
- Golub TR, Slonim DK, Tamayo P, Huard C, Gaasenbeek M, Mesirov JP, Coller H, Loh ML, Downing JR, Caligiuri MA, Bloomfield CD, Lander ES: Molecular classification of cancer: class discovery and class prediction by gene expression monitoring. Science 1999, 286(5439):531–7. 10.1126/science.286.5439.531View ArticlePubMedGoogle Scholar
- Gentleman R: Using GO for statistical analyses. In Compstat 2004 Proceedings in Computational Statistics. Jaromir Antoch: Physica Verlag; 2004:171–80.View ArticleGoogle Scholar
- Yeoh EJ, Ross ME, Shurtleff SA, Williams WK, Patel D, Mahfouz R, Behm FG, Raimondi SC, Relling MV, Patel A, Cheng C, Campana D, Wilkins D, Zhou X, Li J, Liu H, Pui CH, Evans WE, Naeve C, Wong L, Downing JR: Classification, subtype discovery, and prediction of outcome in pediatric acute lymphoblastic leukemia by gene expression profiling. Cancer Cell 2002, 1(2):109–10. 10.1016/S1535-6108(02)00032-6View ArticleGoogle Scholar
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.