On reliable discovery of molecular signatures
© Nilsson et al; licensee BioMed Central Ltd. 2009
Received: 02 April 2008
Accepted: 29 January 2009
Published: 29 January 2009
Molecular signatures are sets of genes, proteins, genetic variants or other variables that can be used as markers for a particular phenotype. Reliable signature discovery methods could yield valuable insight into cell biology and mechanisms of human disease. However, it is currently not clear how to control error rates such as the false discovery rate (FDR) in signature discovery. Moreover, signatures for cancer gene expression have been shown to be unstable, that is, difficult to replicate in independent studies, casting doubts on their reliability.
We demonstrate that with modern prediction methods, signatures that yield accurate predictions may still have a high FDR. Further, we show that even signatures with low FDR may fail to replicate in independent studies due to limited statistical power. Thus, neither stability nor predictive accuracy are relevant when FDR control is the primary goal. We therefore develop a general statistical hypothesis testing framework that for the first time provides FDR control for signature discovery. Our method is demonstrated to be correct in simulation studies. When applied to five cancer data sets, the method was able to discover molecular signatures with 5% FDR in three cases, while two data sets yielded no significant findings.
Our approach enables reliable discovery of molecular signatures from genome-wide data with current sample sizes. The statistical framework developed herein is potentially applicable to a wide range of prediction problems in bioinformatics.
Unfortunately, existing computational approaches often fail to distinguish between the different objectives of prediction and discovery. If molecular signatures are to be used for discovery, then the primary objective is to control the false discovery rate (FDR) with respect to the optimal (true) signature. On the other hand, if the end goal is an accurate predictor, then the FDR of the gene signature is not important in itself. However, it has hitherto not been possible to directly address FDR control, since an operational definition of the optimal signature (a "gold standard") has not been available. Therefore, current methods for signature discovery resort to optimizing prediction accuracy, implicitly assuming that the FDR is thereby kept reasonably low, even though there is no a priori reason to assume that this is the case. Recently, the stability of a signature, that is, the expected overlap between signatures derived from replicated experiments, has been suggested as an alternative quality measure [7, 8]. Signatures derived from cancer gene expression data have been found to be unstable, raising concerns that existing signature discovery methods may not be sound [9, 10]. While the stability measure seems intuitively reasonable and cleverly avoids the gold standard problem, it has not been shown that low stability actually indicates high FDR.
In this paper, we build upon a recently discovered operational definition of the optimal signature to study the actual FDR in signature discovery. First, we demonstrate that high FDR can occur even with very accurate predictors. Therefore, current methods for signature discovery that focus on optimizing prediction accuracy offer no means of controlling the FDR. Second, we show that signatures can be highly unstable even when the FDR is kept low. Thus, reliable signature discovery may be possible in spite of the recent reports of unstable signatures in cancer [9, 10]. Third, we propose a novel hypothesis testing procedure based on our definition of the optimal signature that for the first time directly addresses signature FDR. We show that our method achieves FDR control on simulated data. Application to well-known cancer data sets uncovers three novel molecular signatures for leukemia, colon and breast cancer.
The optimal signature
where gS denotes a predictor on the subspace S of corresponding to the variable set S. Unfortunately, this criterion does not yield a unique S* in general, and there are examples of data distributions such that no tractable (polynomial-time) algorithms exist for computing S* [, pp. 562]. Consequently, most research has focused on heuristic algorithms for discovering approximate signatures with near-optimal prediction accuracy .
That is, S* consists precisely of the variables i such that the error probability of the optimal predictor g* increases when i is removed. The required restriction is that the data density f (x) is everywhere strictly positive. This condition is satisfied by nearly all commonly used statistical models, including the exponential family, and we believe that it is reasonable for biological data. A formal proof of the correctness of (2) is given in Additional File 1.
Note that S* may be quite different from the set of variables that are marginally correlated with the phenotype (e.g., differentially expressed genes). This is because some correlated variables may be "redundant" for prediction: while these do contain information about the phenotype, that information is also present in other variables, so that the redundant variables are excluded from S*. Indeed, it can be proved that S* only contains variables X i that are conditionally dependent on Y regardless of what other variable set is conditioned on . In this sense, S* constitutes the variables "directly" correlated with Y. Moreover, some variables may be predictive only when considered together with certain other variables in a multivariate fashion, and thus S* may contain variables that are not detectable by standard univariate methods .
The simple form (2) immediately suggests a general, linear-time, asymptotically correct algorithm for discovering S* from data, as described elsewhere . Here, we make use of the fact that (2) permits S* to be calculated for any given data distribution, thus providing the gold standard required for evaluating signature discovery methods and developing hypothesis testing procedures.
Accurate predictions despite high signature FDR
The likely reason for this behavior is that modern predictive methods such as the SVM have internal mechanisms for suppressing noise (regularization). They are therefore rather insensitive to false positives within the signature. For prediction purposes, it is more important that the signature does contain some true positives genes, while a large fraction of irrelevant genes may be tolerated without degrading predictive accuracy. As a consequence, discovering signatures by optimizing prediction accuracy should not be expected control FDR, as we will further demonstrate below.
Unstable Signatures with Low FDR
A Statistical Framework for Signature Discovery
where denote the weights of the optimal classifier. Assuming that the classifier used is consistent, we have that [w i ] → as sample size increases. Hence, in this case we can equivalently test the null hypothesis [w i ] = 0. More complicated parametric forms such as polynomials in x i could be used in a similar way, although the number of weights would increase accordingly.
Since the statistical distribution of w i is unknown, we used a bootstrap technique to test whether [w i ] = 0. By sampling with replacement from the given data set and re-training the classifier on each sample, we obtain B vectors . For each variable i, the corresponding are then used to obtain a bootstrap confidence interval for w i . This interval is inverted to obtain a bootstrap p-values p i for each variable i (that is, the null hypothesis is rejected at level α if the (1 - α) confidence interval does not cover zero). Importantly, this procedure preserves the full dependency structure of the data distribution. Finally, FDR control was performed using the Benjamini-Hochberg procedure .
We repeated the simulation study using the popular Recursive Feature Elimination (RFE) method  to discover signatures. While this method did produce accurate predictive models (data not shown), we observed that FDR was high (above 40% in this experiment) and depended on the effect size in an unpredictable manner. Indeed, optimizing prediction accuracy by RFE does not guarantee a reliable signature. High FDR was also observed when choosing the signature S as a fixed-size "top list" by the rank according to the w g statistics (Figure 4D). We have also previously observed high FDR for other methods that optimize the signature for prediction accuracy . Often, these methods attempt to include more variables in the signature when the prediction problem is harder, thus sacrificing FDR control for better predictive accuracy. Conversely, for less difficult prediction problems, many true positives may be removed from the signature because they do not influence predictive power discernably.
Application to Cancer Gene Expression
Results on cancer gene expression data
Data set (ref.)
97.0 ± 4.2
92.6 ± 3.0
81 ± 7.2
65 ± 4.3
van't Veer (5)
62 ± 8.4
As a negative control, we applied our bootstrap test on randomized versions of each original data set where the phenotype values were randomly permuted, corresponding to the complete null hypothesis. This yielded zero significant genes in each case, confirming that we do not obtain spurious findings. In contrast, when applying the RFE method to randomized data, we consistently obtained even larger signatures than with the real data sets. We also tested each signature on an independent data set, confirming that the signatures are indeed predictive.
For comparison, we performed a conventional differential expression test for each data set using the t-test statistic with the Benjamini-Hochberg correction (Table 1). This identified a substantially larger set of genes than the bootstrap method – in one case, more than half of the genes tested were significant. This illustrates the ability of the gene signature approach to distinguish the genes directly related to the phenotype variable from a much larger set of differentially expressed genes: many of the latter turn out to be "redundant" for prediction, meaning that they are correlated with the phenotype only indirectly, through genes in S*.
Molecular signatures offer a systematic way to focus on the genes most directly associated with a given phenotype, and may yield valuable insights into the underlying biological system. It is therefore unfortunate that the reliability of signatures per se is poorly understood. Since no gold standard for signature discovery has been available, validation of discovered signatures often amounts to mining the scientific literature for documented connections between the phenotype being studied and the elements (genes) of a hypothesized signature. However, this approach is necessarily biased and rather speculative: it is by no means clear that a gene should be included in a predictive signature simply because it is somehow "related" to the phenotype. For example, approximately 25% of all known human genes have some documented relation to cancer , but it is unlikely that all of these should be included in an optimal signature for cancer prediction.
To address this issue, we have herein developed a statistical method for signature discovery based on a formal definition of the "gold standard" optimal signature. This allows for assessing the reliability of signatures without detailed knowledge of the biological system. To our knowledge, our method is the first to provide statistical guarantees for the reliability of molecular signatures, although we note that random forests are similar to our bootstrap testing scheme and also give indications of what variables are important for prediction.
For two of the gene expression data sets investigated, including the well-studied cancer data by van't Veer et al. , our method did not call any genes significant, indicating that these data sets did not contain sufficient information to uncover gene signatures at the specified false discovery rate (5%). We emphasize that this does not necessarily mean that it is infeasible to construct predictive models for these studies, but merely that it is difficult to determine which genes are responsible for the predictive accuracy. In this sense, discovering reliable gene signatures can be a harder problem than obtaining accurate predictors. Prediction and signature discovery are two separate problems, and must be treated differently.
For simplicity, we have here restricted our analysis to two-class problems and linear predictors. However, the proposed method is applicable to any learning method for which a reasonably well-powered statistic can be derived to test the signature null hypothesis. Continuous phenotype variables can easily be addressed by substituting the classification methods used herein for regression methods, such as ridge regression  or the relevance vector machine . General methods for handling non-linear dependencies have also been described [13, 24], although it is unclear whether signature discovery from gene expression data would benefit from these more complex models with currently available sample sizes.
Some technical issues remain to be considered. First, testing the null hypothesis [w i ] = 0 is technically correct only in the limit of large samples where [w i ] → . While our simulation studies indicate correct behavior for the sample sizes tested, this issue warrants further study. Second, bootstrap hypothesis testing is known to provide only approximate p-values, satisfying the inequality P(p ≤ α) ≤ α + (1/l), where l is the sample size . While the additional term (1/l) was negligible in our simulations, this should be verified in each particular case before applying bootstrap testing. A possible future improvement could be to estimate this term from simulations and correct the bootstrap p-values accordingly, thereby "calibrating" the method.
As we have shown, neither predictive accuracy nor stability constitute valid measures of FDR for molecular signatures. Indeed, highly accurate predictions may be obtained despite an FDR as high as 50% (Figure 2), while in situations where many weak effects are present and statistical power is low, signatures can be unstable at an FDR as low as 2.5% (Figure 3). This analysis explains at least some of the difficulties with reproducing cancer gene expression signatures [7, 8] and possibly also the similar reproducibility problems of recent association studies in complex diseases . Moreover, it suggests that this lack of reproducibility need not be problematic.
We have developed and validated a statistical hypothesis testing framework that for the first time provides false discovery rates control for signature discovery. In application to cancer gene expression, we have showed that reliable signature discovery is feasible with currently available sample sizes. Many important problems in bioinformatics are prediction problems and may benefit from reliable signature discovery. We therefore hope that our method will be of general interest.
where denotes statistical expectation. The stability is always between 0 (no expected overlap) and 1 (complete overlap).
where is the error function.
To evaluate signature error rates, we used the fact that for f (x | y) = N (yμ, Σ), the optimal separating hyperplane has normal vector w* = Σ-1μ, and so the optimal set S* can be determined as the nonzero components of this vector.
For hypothesis testing, we used a parametric bootstrap with B = 50 repetitions, fitting a Gaussian distribution N (μ i , σ i ) to the observed prior to computing two-sided p-values. In preliminary studies, the difference between this method and a nonparametric bootstrap with B = 1000 was negligible, while the parametric version is computationally more efficient since a much smaller B can be used. The SVM , KFD  and VW  methods were implemented as previously described. In all experiments, the SVM C-parameter and the KFD regularization parameter were set to 1. Recursive Feature Elimination (RFE) was performed as previously described , using the radius-margin bound  as accuracy measure and removing 20% of the genes in each iteration.
Microarray data sets [1–5] were preprocessed by removing genes displaying small variation, keeping the 5,000 most variable genes in each case, except for the data sets by van't Veer et al.  and Alon et al.  which were preprocessed in a similar fashion by the original authors. Genes were normalized to zero mean and unit standard deviation prior to SVM training, following standard practise for kernel methods. Independent test data sets [27–29] were normalized in the same fashion. No other preprocessing was done prior to classifier training or testing.
where Acc+ and Acc- are the accuracy measures for each class. Except for the independent test sets, these were measured by cross-validation, where in each round a randomized set consisting of 2/3 of the samples was used for training, and the remaining 1/3 was used for testing. Splits were balanced so that class frequencies were equal between training/test data. Mean and standard deviation of the balanced test error over 50 cross-validation repetitions are reported.
The authors would like to thank Drs. José M. Peña and Albert Compte for helpful discussions. This work was supported by grants from the Ph.D. Programme in Medical Bioinformatics at Karolinska Institutet (RN), Clinical Gene Networks AB, Vinnova (JT), Swedish Research Council (JT) and Linköping University.
- Alon U, Barkai N, Notterman DA, Gish K, Ybarra S, Mack D, Levine AJ: Broad patterns of gene expression revealed by clustering analysis of tumor and normal colon tissues probed by oligonucleotide arrays. Proc Natl Acad Sci U S A 1999, 96(12):6745–6750. 10.1073/pnas.96.12.6745PubMed CentralView ArticlePubMedGoogle Scholar
- Golub TR, Slonim DK, Tamayo P, Huard C, Gaasenbeek M, Mesirov J, Coller H, Loh M, Downing J, Caligiuri M, Bloomfield C, Lander ES: Molecular classifiation of cancer: class discovery and class prediction by gene expression monitoring. Science 1999, 286(5439):531–537. 10.1126/science.286.5439.531View ArticlePubMedGoogle Scholar
- Singh D, Febbo PG, Ross K, Jackson DG, Manola J, Ladd C, Tamayo P, Renshaw AA, D'Amico AV, Richie JP, Lander ES, Loda M, Kanto3 PW, Golub TR, Sellers WR: Gene expression correlates of clinical prostate cancer behavior. Cancer Cell 2002, 1(2):203–209. 10.1016/S1535-6108(02)00030-2View ArticlePubMedGoogle Scholar
- van't Veer LJ, Dai H, Vijver MJ, He YD, Hart AAM, Mao M, Peterse HL, Kooy K, Marton MJ, Witteveen AT, Schreiber GJ, Kerkhoven RM, Roberts C, Linsley PS, Bernards R, Friend SH: Gene expression profiling predicts clinical outcome of breast cancer. Nature 2002, 415(6871):530–536. 10.1038/415530aView ArticleGoogle Scholar
- Wang Y, Klijn JG, Zhang Y, Sieuwerts AM, Look MP, Yang F, Talantov D, Timmermans M, Meijer-van Gelder ME, Yu J, Jatkoe T, Berns EM, Atkins D, Foekens JA: Gene-expression profiles to predict distant metastasis of lymph-node-negative primary breast cancer. Lancet 2005, 365(9460):671–679.View ArticlePubMedGoogle Scholar
- Bogaerts J, Cardoso F, Buyse M, Braga S, Loi S, Harrison JA, Bines J, Mook S, Decker N, Ravdin P, Therasse P, Rutgers E, van't Veer LJ, Piccart M: Gene signature evaluation as a prognostic tool: challenges in the design of the MINDACT trial. Nat Clin Pract Oncol 2006, 3(10):540–551. 10.1038/ncponc0591View ArticlePubMedGoogle Scholar
- Michiels S, Koscielny S, Hill C: Prediction of cancer outcome with microarrays: a multiple random validation strategy. Lancet 2005, 365(9458):488–492. 10.1016/S0140-6736(05)17866-0View ArticlePubMedGoogle Scholar
- Ein-Dor L, Kela I, Getz G, Givol D, Domany E: Outcome signature genes in breast cancer: is there a unique set? Bioinformatics 2005, 21(2):171–178. 10.1093/bioinformatics/bth469View ArticlePubMedGoogle Scholar
- Ein-Dor L, Zuk O, Domany E: Thousands of samples are needed to generate a robust gene list for predicting outcome in cancer. Proc Natl Acad Sci U S A 2006, 103(15):5923–5928. 10.1073/pnas.0601231103PubMed CentralView ArticlePubMedGoogle Scholar
- Ioannidis JP: Microarrays and molecular research: noise discovery? Lancet 2005, 365(9458):454–455.View ArticlePubMedGoogle Scholar
- Devroye L, Györfi L, Lugosi G: A probabilistic theory of pattern recognition. Applications of mathematics. New York: Springer-Verlag; 1996.View ArticleGoogle Scholar
- Guyon I, Elisseeff A: An introduction to variable and feature selection. Journ Mach Learn Res 2003, 3: 1157–1182. 10.1162/153244303322753616Google Scholar
- Nilsson R, Peña JM, Björkegren J, Tegnér J: Consistent feature selection for pattern recognition in polyomial time. Jour of Mach Learn Res 2007, 8: 589–612.Google Scholar
- Nilsson R, Peña JM, Björkegren J, Tegnér J: Detecting multivariate differentially expressed genes. BMC Bioinformatics 2007, 8: 150. 10.1186/1471-2105-8-150PubMed CentralView ArticlePubMedGoogle Scholar
- Benjamini Y, Hochberg Y: Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Statist Soc B 1995, 57: 289–300.Google Scholar
- Frayling TM: Genome-wide association studies provide new insights into type 2 diabetes aetiology. Nat Rev Genet 2007, 8(9):657–662. 10.1038/nrg2178View ArticlePubMedGoogle Scholar
- Schäfer J, Strimmer K: An empirical Bayes approach to inferring large-scale gene association networks. Bioinformatics 2005, 21(6):754–764. 10.1093/bioinformatics/bti062View ArticlePubMedGoogle Scholar
- Cortes C, Vapnik V: Support-Vector Networks. Mach Learn 1995, 20(3):273–297.Google Scholar
- Mika S, Ratsch G, Weston J, Scholkopf B, Muller K: Fisher Discriminant Analysis with Kernels. In Proceedings of IEEE Neural Networks for Signal Processing Workshop Edited by: Hen YH, Larsen J, Wilson E. 1999, 41–48.Google Scholar
- Guyon I, Weston J, Barnhill S, Vapnik V: Gene Selection for Cancer Classification using Support Vector Machines. Mach Learn 2002, 46: 389–422. 10.1023/A:1012487302797View ArticleGoogle Scholar
- Nilsson R, Peña JM, Björkegren J, Tegnér J: Evaluating feature selection for SVMs in high dimensions. Proceedings of the 17th European Conference on Machine Learning 2006, 719–726.Google Scholar
- Heorl A, Kennard R: Ridge regression: biased estimation of nonorthogonal problems. Technometrics 1970, 12: 69–82. 10.2307/1267352View ArticleGoogle Scholar
- Tipping ME: Sparse Bayesian learning and the relevance vector machine. Journ Mach Learn Res 2001, 1: 211–244. 10.1162/15324430152748236Google Scholar
- Li F, Yang Y, Xing EP: From LASSO regression to feature vector machine. In Advances in Neural Information Processing Systems 18. Edited by: Weiss Y. MIT Press, Cambridge; 2005:411–418.Google Scholar
- Efron B, Tibshirani RJ: An introduction to the bootstrap. Chapman & Hall, Inc. New York; 1993.View ArticleGoogle Scholar
- Vapnik VN: Statistical Learning Theory. John Wiley and Sons, Inc. New Jersey; 1998.Google Scholar
- Yu Y, Landsittel D, Jing L, Nelson J, Ren B, Liu L, McDonald C, Thomas R, Dhir R, Finkelstein S, Michalopoulos G, Becich M, Luo JH: Gene expression alterations in prostate cancer predicting tumor aggression and preceding development of malignancy. J Clin Oncol 2004, 22(14):2790–2799. 10.1200/JCO.2004.05.158View ArticlePubMedGoogle Scholar
- Campo Dell'Orto M, Zangrando A, Trentin L, Li R, Liu W, te Kronnie G, Basso G, Kohlmann A: New data on robustness of gene expression signatures in leukemia: comparison of three distinct total RNA preparation procedures. BMC Genomics 2007, 8: 188. 10.1186/1471-2164-8-188PubMed CentralView ArticlePubMedGoogle Scholar
- Available from the NCBI Gene Expression Omnibus, accession GSE10960[http://www.ncbi.nlm.nih.gov/geo/]
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.