Sparse non-negative matrix factorizations via alternating non-negativity-constrained least squares for microarray data analysis (original) (raw)

Journal Article

,

College of Computing, Georgia Institute of Technology, 266 Ferst Drive, Atlanta, GA 30332, USA

*To whom correspondence should be addressed.

Search for other works by this author on:

College of Computing, Georgia Institute of Technology, 266 Ferst Drive, Atlanta, GA 30332, USA

*To whom correspondence should be addressed.

Search for other works by this author on:

Received:

02 November 2006

Revision received:

19 February 2007

• PDF
• Split View
• • Article contents
  • Figures & tables
  • Video
  • Audio
  • Supplementary Data
• Cite

Cite

Hyunsoo Kim, Haesun Park, Sparse non-negative matrix factorizations via alternating non-negativity-constrained least squares for microarray data analysis, Bioinformatics, Volume 23, Issue 12, June 2007, Pages 1495–1502, https://doi.org/10.1093/bioinformatics/btm134
Close

Navbar Search Filter Mobile Enter search term Search

Abstract

Motivation: Many practical pattern recognition problems require non-negativity constraints. For example, pixels in digital images and chemical concentrations in bioinformatics are non-negative. Sparse non-negative matrix factorizations (NMFs) are useful when the degree of sparseness in the non-negative basis matrix or the non-negative coefficient matrix in an NMF needs to be controlled in approximating high-dimensional data in a lower dimensional space.

Results: In this article, we introduce a novel formulation of sparse NMF and show how the new formulation leads to a convergent sparse NMF algorithm via alternating non-negativity-constrained least squares. We apply our sparse NMF algorithm to cancer-class discovery and gene expression data analysis and offer biological analysis of the results obtained. Our experimental results illustrate that the proposed sparse NMF algorithm often achieves better clustering performance with shorter computing time compared to other existing NMF algorithms.

Availability: The software is available as supplementary material.

Contact: hskim@cc.gatech.edu, hpark@acc.gatech.edu

Supplementary information: Supplementary data are available at Bioinformatics online.

1 INTRODUCTION

In many data-mining problems, dimension reduction is imperative for efficient manipulation of massive quantity of high-dimensional data. The subspace method has demonstrated its success in numerous pattern recognition tasks including efficient classification (Kim et al., 2005), clustering (Ding et al., 2002) and fast search (Berry et al., 1999). There are two general approaches for reducing dimensionality, i.e. feature extraction and feature selection. Feature extraction is transforming the existing features into a lower dimensional space, while feature selection is selecting a subset of the existing features without a transformation. For feature extraction, principal component analysis (PCA), linear discriminant analysis (LDA) and non-negative matrix factorization (NMF) have been widely used. Many practical pattern recognition problems require non-negativity constraints. For example, pixels in digital images and chemical concentrations in bioinformatics are non-negative. NMF is a useful technique in approximating these high-dimensional data.

Given a non-negative matrix A of size ⁠, where each column of A corresponds to a data point in the _m_-dimensional space, and a positive integer ⁠, NMF finds two non-negative matrices and so that

(1)

A solution to the NMF problem can be obtained by solving the following optimization problem:

(2)

where is a basis matrix, is a coefficient matrix, is the Frobenius norm and W, H ≥ 0 means that all elements of W and H are non-negative. Due to ⁠, dimension reduction is achieved and a lower dimensional representation of A in a _k_-dimensional space is given by H.

Since NMF may give us direct interpretation due to non-subtractive combinations of non-negative basis vectors, it has recently received much attention and it has been applied to many interesting problems including text data mining (Chagoyen et al., 2006; Lee and Seung, 1999; Pauca et al., 2004) gene expression data analysis (Brunet et al., 2004; Carmona-Saez et al., 2006; Gao and Church, 2005; Kim and Tidor, 2003; Maher et al., 2006), microarray comparative genomic hybridization (aCGH) data (Carrasco et al., 2006) and functional characterization of gene lists (Pehkonen et al., 2005). One of the interesting properties of NMF is that it often generates sparse basis vectors that allow us to discover parts-based basis vectors. NMF generated holistic basis images instead of parts-based basis images for a facial image dataset in the results presented in Li et al. (2001) and Hoyer (2004). Several approaches (Dueck et al., 2005; Hoyer, 2004; Pascual-Montano et al., 2006; Pauca et al., 2006) have been proposed to explicitly control the degree of sparseness of W and H. However, if strong sparsity constraints are imposed on the basis matrix W, some useful information for gene selection in microarray data analysis may be lost. Some other approaches (Gao and Church, 2005; Pauca et al., 2004) imposed sparsity constraints on only the coefficient matrix H, but there may be mathematical difficulties related to convergence.

In this article, we introduce a sparse NMF algorithm that can control the degree of sparseness in the coefficient matrix H via alternating non-negativity-constrained least squares, and apply it to microarray data analysis. The rest of this article is organized as follows. We give brief overviews on various NMF algorithms in Section 2.1. In Section 2.2, we introduce sparse NMF formulations and algorithms based on alternating non-negativity-constrained least squares and discuss their convergence properties. In Section 3, we describe stopping criteria, and present experimental results and biological analysis illustrating properties of the proposed sparse NMF. Summary is given in Section 4.

2 METHODS

2.1 Review of NMF algorithms

Lee and Seung (1999, 2000) suggested NMF algorithms based on multiplicative update rules of W and H. The distance is non-increasing under the following multiplicative update rules:

for 1 ≤ q ≤ k and 1 ≤ j ≤ n,

for 1 ≤ i ≤ m and 1 ≤ q ≤ k. The divergence D(A, WH) = is non-increasing under the different multiplicative update rules (Lee and Seung, 2000). Gonzales and Zhang (2005) claimed that these non-increasing properties of multiplicative update rules may not imply the convergence to a stationary point within realistic amount of runtime for problems of meaningful sizes. More detailed review on NMF algorithms can be found in Berry et al. (2006).

Hoyer (2004) devised a sparse NMF algorithm based on the projected gradient descent method (SNMF/PGD) in order to constrain NMF to find solution with desired sparseness in W and H. To impose sparseness constraints on only one matrix W or H, this algorithm uses a multiplicative update rule for updating the counter matrix, which suffers from slow convergence. Puscual-Montano et al. (2006) claimed that non-smooth NMF (_ns_NMF) outperformed previous sparse NMF variants for their synthetic and real datasets. The _ns_NMF (Puscual-Montano et al., 2006) is also based on multiplicative update rules.

Pauca et al. (2006) proposed a constrained NMF (CNMF) formulation,

(3)

where α and β are regularization parameters. A sparse NMF algorithm using the following least squares,

(4)

has appeared in Pauca et al. (2004) and Gao and Church (2005). This algorithm sets negative values in H to zero during iterations. However, setting negative values to zero for imposing non-negativity is not recommended, since there is no guarantee that the algorithm converges (Bro and de jong, 1997). The parameter β in Equation (4) has a scaling effect since a large value of β would suppress ⁠. As is suppressed, may grow relatively large, and therefore, the algorithm needs column normalization of W during iterations. However, the normalization of W changes the objective function, and this makes convergence analysis difficult. It is well known that a quadratic penalty corresponds to Gaussian priors and does not encourage sparsity but rather scales the result giving non-sparse low values. Thus, _L_1-norm based formulations would be more appropriate than _L_2-norm based formulations so as to control sparsity (Tibshirani, 1996).

2.2 Sparse NMFs based on alternating non-negativity-constrained least squares

In order to enforce sparseness on W or H in the NMF presented in Equation (1), we introduce two formulations and the corresponding algorithms for sparse NMFs, i.e. SNMF/L for sparse W (where ‘L’ denotes the sparseness imposed on the left factor) and SNMF/R for sparse H (where ‘R’ denotes the sparseness imposed on the right factor). Our sparse NMF formulations that impose the sparsity on a factor of NMF utilize _L_1-norm minimization and the corresponding algorithms are based on alternating non-negativity constrained least squares (ANLS). Each sub-problem is solved by a fast non-negativity constrained least squares (NLS) algorithm (van Benthem and Keenan, 2004) that is improved upon the active set based NLS method. Bro and de Jong (1997) made a substantial speed improvement to Lawson and Hanson's algorithm (Lawson and Hqnson, 1974) for multiple right-hand-side cases. van Benthem and Keenan (2004) devised an algorithm that further improves the performance of NLS.

2.2.1 Formulations for Sparse NMFs

SNMF/R: To apply sparseness constraints on H, we formulate the following SNMF/R optimization problem:

(5)

where is the _j_-th column vector of H, is a parameter to suppress ⁠, and is a regularization parameter to balance the trade-off between the accuracy of the approximation and the sparseness of H. The SNMF/R algorithm begins with the initialization of W with non-negative values. Then, it iterates the following ANLS until convergence:

(6)

where is a row vector with all components equal to one and is a zero vector, and

(7)

where Ik is an identity matrix of size and 0_k_ × m is a zero matrix of size ⁠. Equation (6) minimizes _L_1-norm of columns of which imposes sparsity on H.

SNMF/L: To impose sparseness constraints on W, we introduce the SNMF/L formulation:

(8)

where W (i, :) is the _i_-th row vector of W, is a parameter to suppress ⁠, and is a regularization parameter to balance the trade-off between the accuracy of the approximation and the sparseness of W. The algorithm for SNMF/L is also based on ANLS.

2.2.2 Convergence properties of sparse NMF algorithms

We show the convergence property of the sparse NMF algorithms. Since the convergence properties of SNMF/L and SNMF/R are essentially the same, we will only discuss the case of SNMF/R in more detail. Under conditions of ⁠, ⁠, and due to ⁠, Equation (5) can rewritten as

(9)

which is differentiable with respect to W or H. SNMF/R contains two sub-problems for two-block minimization scheme. Grippo and Sciandrone (2000) showed that the two-block coordinate descent method does not require each sub-problem to have a unique solution for convergence. The objective function in Equation (9) is coercive on the feasible set; as the feasible set is closed, the intersection of any level set of this function with the feasible set is compact. Therefore, any minimization process that reduces the objective function and preserves feasibility generates points that remain in a compact set. The existence of accumulation points and the differentiability of the objective function in Equation (9) imply that the assumptions of Grippo and Sciandrone's Corollary (Grippo and Sciandrone, 2000) are satisfied, so that we can establish that the two-block minimization process is convergent, in the sense that every accumulation point is a critical point of the problem shown in Equation (9). Similarly, it can be shown that the algorithm SNMF/L converges to a stationary point.

3 RESULTS

3.1 Datasets description

We used the leukemia gene expression dataset (ALLAML) (Golub et al., 1999) and the central nervous system tumors dataset (CNS) (Pomeroy et al., 2002). The ALLAML dataset contains acute lymphoblastic leukemia (ALL) that has B and T cell subtypes, and acute myelogenous leukemia (AML) that occurs more commonly in adults than in children. This gene expression dataset consists of 38 bone marrow samples (19 ALL-B, 8 ALL-T and 11 AML) with 5000 genes. The central nervous system dataset is composed of four categories of CNS tumors with 5597 genes. It consists of 34 samples representing four distinct morphologies: 10 classic medulloblastomas, 10 malignant gliomas, 10 rhabdoids and 4 normal cerebella. All datasets used by us contain only non-negative entries. All algorithms were implemented in Matlab 6.5 (MATLAB, 1992). The Matlab codes for NMF using divergence-based multiplicative update rules were obtained from Brunet et al. (2004) and modified to implement _ns_NMF (Puscual-Montano et al., 2006). All our experiments were performed on a P3 600 MHz machine with 512 MB memory.

3.2 Biclustering

We applied the non-negative factorization of Equation (1) to perform clustering analysis of a data matrix. The rows of a microarray data matrix A represent genes and the columns experiments. We can use the basis matrix W to divide the m genes into k gene-clusters and the coefficient matrix H to divide the n samples into k sample-clusters. Typically, gene i is assigned to gene-cluster q if the W(i, q) is the largest element in W (i, :) and sample j is assigned to sample-cluster q if the H(q, j) is the largest element in H(:, j).

3.3 Stopping criteria

NMF using divergence-based multiplicative update rules (i.e. NMF/DUR (Brunet et al., 2004) and _ns_NMF (Pascual-Montano et al., 2006)) in our implementation stops if has not changed for more than 40 convergence tests (each made at 10 iterations), where is the connectivity matrix of which entry is if samples i and j belong to the same sample-cluster, and if they belong to different sample-clusters.

For SNMF/L and SNMF/R, we tested convergence at every five iterations by the combined convergence criterion using the Karush-Kuhn-Tucker (KKT) optimality conditions and the convergence of positions of the largest elements in rows of W and columns of H. Our sparse NMF algorithms stop if both the positions of the largest elements in rows of W, i.e. ⁠, and the positions of the largest elements in columns of H, i.e. ⁠, have not changed for more than or equal to 10 convergence tests, where is the position of the largest element in the _i_-th row of W and is the positions of the largest element in the _j_-th column of H, and the following KKT conditions are satisfied. The KKT conditions for each objective function f(W, H) with non-negativity constraints W ≥ 0 and H ≥ 0 are

(10)

These conditions can be rewritten as

(11)

where the minimum is taken component wise (Gonzales and Zhang, 2005). Let Δ o be the KKT residual measured by the _L_1-vector norm,

(12)

We count the number of the elements in W that did not converge yet, i.e. δ W = #(min(W, /∂W) ≠ 0), and the number of the elements in H that did not converge yet, i.e. δ H = #(min(H, /∂H) ≠ 0). We define the following normalized KKT residual:

(13)

which reflects the average of convergence errors for elements in W and H that did not converge. The mathematical convergence criterion is defined as

(14)

where Δ1 is the Δ value in the first iteration and ɛ is a tolerance. We used for our experiments.

3.4 Clustering performance comparison

To measure the performance of NMFs in clustering, we used purity and entropy. Suppose we are given l categories (true class labels), while NMF generates k clusters. Purity is given by

where n is the total number of samples and is the number of samples in the cluster q that belong to original class j (1 ≤ j ≤ l). The larger the value of purity, the better the clustering performance. Entropy is defined as follows:

where l denotes the number of original class labels and nq is the size of cluster q. The smaller the value of entropy, the better the clustering quality.

The parameter η in Equation (5) is important in keeping small. We set it to be the square of the maximal element in A. For the initialization of W in SNMF/R, the elements in the initial matrix W were randomly chosen and normalized so that the columns of the basis matrix W have unit _L_2-norm, i.e. ⁠.

Tables 1 and 2 show the results of SNMF/R with various values of β on the ALLAML dataset with k = 3 and on the CNS tumors dataset with k = 4, respectively. We compared our proposed SNMF/R with NMF based on divergence-based multiplicative update rules (Brunet et al., 2004; Lee and Seung, 2000). Average sparseness, purity and entropy were computed by running each algorithm five times with different random initializations. By increasing β in SNMF/R, we could obtain a sparser H. SNMF/R algorithm achieved better clustering performance (higher purity, lower entropy) than NMF/DUR within certain range of β with shorter computing time. By increasing α in SNMF/L, we could enhance the sparsity of W (results are not shown). SNMF/L can be applied to obtain parts-based basis vectors.

Table 1.

Performance dependency of SNMF/R (k=3) on various β values on the leukemia data matrix of size 5000 × 38. We present the average percentages of zero elements in W and H over five runs with different random initializations. We also present average purity, entropy, computing time (in seconds) and the number of iterations

*For NMF using divergence-based multiplicative update rules (NMF/DUR) \citep {brunet04metagenes, the average percentages of the number of very small non-negative elements that are smaller than 10−8 in W and H are presented.

Table 1.

Performance dependency of SNMF/R (k=3) on various β values on the leukemia data matrix of size 5000 × 38. We present the average percentages of zero elements in W and H over five runs with different random initializations. We also present average purity, entropy, computing time (in seconds) and the number of iterations

*For NMF using divergence-based multiplicative update rules (NMF/DUR) \citep {brunet04metagenes, the average percentages of the number of very small non-negative elements that are smaller than 10−8 in W and H are presented.

Table 2.

Performance dependency of SNMF/R (k=4) on various β values on the CNS tumors data matrix of size 5597 × 34. We present the average percentages of zero elements in W and H over five runs with different random initializations. We also present average purity, entropy, computing time (in seconds) and the number of iterations

*For NMF using divergence-based multiplicative update rules (NMF/DUR) (Brunet et al., 2004), the average percentages of the number of very small non-negative elements that are smaller than 10−8 in W and H are presented.

Table 2.

Performance dependency of SNMF/R (k=4) on various β values on the CNS tumors data matrix of size 5597 × 34. We present the average percentages of zero elements in W and H over five runs with different random initializations. We also present average purity, entropy, computing time (in seconds) and the number of iterations

*For NMF using divergence-based multiplicative update rules (NMF/DUR) (Brunet et al., 2004), the average percentages of the number of very small non-negative elements that are smaller than 10−8 in W and H are presented.

We compared our methods with other sparse NMF variants on the ALLAML dataset. We tested Hoyer's sparse NMF based on the projected gradient descent method by his Matlab implementation with the sparseness control parameter to impose sparsity constraints on only H (see Hoyer (2004) for the details). The average percentages of zero elements in W and H obtained from SNMF/PGD were 0.12 and 21.75%, respectively. SNMF/R showed significantly better clustering performance than SNMF/PGD (average purity = 0.895 and average entropy = 0.280). Moreover, SNMF/R required much shorter computing time and smaller number of iterations than SNMF/PGD (average computing time = 110.1 and average iteration number = 517.8). We tested _ns_NMF with the smoothness control parameter θ = 0.5 suggested in Carmona-Saez et al. (2006) for biclustering of gene expression data due to reasonable results from numerous empirical tests. _ns_NMF generated average purity = 0.963 and average entropy = 0.108. The average percentages of elements in the range of in W and H obtained from _ns_NMF were 8.94 and 25.26%, respectively. It took 79.1 s with 698 iterations on the average. It produced the sparsest W among NMF algorithms we compared. The average percentages of elements in the range of in W and H obtained from _ns_NMF were 65.94 and 29.30%, respectively, which are much higher than 0.12 and 0.35% obtained from NMF/DUR. _ns_NMF enhanced the sparseness of W as well as H simultaneously. However, the additional sparsity constraints on W is not always helpful. For example, some genes are over-expressed in samples that belong to more than one clusters. NMFs typically generate W whose rows corresponding to these genes have large values in more than one factors. One can obtain a sparser W by imposing strong sparsity constraints on W. However, the sparser W loses the original information and may give us wrong information that these genes are over-expressed in only one factor. Consequently, one may inappropriately select these genes from the sparse W since they are considered as factor-specific genes. We also tested the probabilistic sparse matrix factorization (PSMF) (Dueck et al., 2005). PSMF also produced a sparser W. Moreover, our SNMF/R showed better performance than PSMF in terms of purity, entropy and computing speed.

We used the model selection method proposed by Brunet et al. (2004) to determine the number of factors. We ran NMF algorithms 30 times to obtain the average connectivity matrix (i.e. consensus matrix) whose entries reflect the probability that samples i and j belong to the same cluster. To measure the dispersion of a consensus matrix C, we defined the dispersion coefficient (ρ) as

(15)

The value of coefficient is ρ = 1 for a perfect consensus matrix (all entries = 0 or 1) and 0 ≤ ρ < 1 for a scattered consensus matrix. After obtaining ρ k values for various k, we can determine the number of clusters from the maximal ρ k. Figure 1 illustrates the determination of the number of clusters in the CNS tumors dataset. SNMF/R with β = 0.01 found perfect consensus matrices for ⁠. In other words, SNMF/R with β = 0.01 generated H matrices that have the same cluster structure in spite of different random initializations of W for 2 ≤ k ≤ 4. SNMF/R yielded finer consensus matrices (higher ρ k) than NMF/DUR for various k values.

Fig. 1.

CNS tumors clustering by SNMF/R with β = 0.01. (Left) The reordered consensus matrices on the CNS tumors dataset. (Right) The corresponding dispersion coefficients. The dispersion coefficient drops when k increases from 4 to 5, indicating a four-cluster split of the data is more stable than a five-cluster split.

3.5 Biological analysis

Figure 2 presents the matrices W and H obtained from SNMF/R with β = 0.01 on the ALLAML leukemia dataset, which produced the lowest approximation error after five runs with different random initializations of W. A column vector of the coefficient matrix H has the contributions of k biological processes to the gene expression of a sample. From the matrix H, we can recognize that ALL-B is dominated by the first biological process. ALL-T is almost controlled by the second biological process. The third biological process is the major component for AML cluster. The cluster of each sample was determined by the positions of the largest elements in columns of H. Only one sample (the 29th sample, AML_13) was incorrectly assigned to ALL-B.

Fig. 2.

W (basis matrix) and H (coefficient matrix) obtained from SNMF/R with β = 0.01 for the ALLAML leukemia dataset (38 samples: 19 ALL-B, 8 ALL-T and 11 AML) with the 5000 most highly varying genes.

A row vector of the basis matrix W indicates the contributions of a gene to the k biological pathways or processes (i.e. k columns of W). Genes can participate in more than one biological process. It is beneficial to investigate genes that have relatively large coefficient in each biological process. We selected factor-specific genes via the non-negative basis matrix obtained from SNMF/R. We define gene_score for the _i_th gene as

(16)

where is a probability that the _i_-th gene contributes to cluster Ω, i.e. ⁠. The gene_score is a real value within the range of [0, 1]. The higher the gene_score value, the more factor-specific the corresponding gene. By using the gene_scores obtained from W, we ranked genes and chose genes whose gene_scores were higher than ⁠, where and are the median and the median absolute deviation (MAD) of gene_scores respectively, and the maximal values in the corresponding rows of W were larger than the median of all elements in W. Total of 730 genes were selected. The numbers of genes chosen for three factors were 423 genes for ALL-B cluster, 104 genes for ALL-T cluster and 203 genes for AML cluster.

Some chosen genes dominantly contribute to only single biological pathway or process. For instance, MB-1 gene (U05259) is most active in the first process. Transcription factor 7 (T-cell specific, HMG-box) (TCF7, X59871) is active in the second process, which is also known as T cell factor-1 (TCF-1). Some genes play a major role in the third process, for example, Interleukin 8 (IL8, M28130), DF D component of complement (adipsin) (CFD, M84526), Cystatin C (amyloid angiopathy and cerebral hemorrhage (CST3, M27891), Chemokine (C-X-C motif) ligand 2 (CXCL2, M57731), etc. Chemokine is a type of cytokines that bind to a specific cell-surface receptor and is critical to the functioning of both innate and adaptive immune responses. Total of 37 genes including MB-1, IL8, CFD and CST3 were the same genes as those found in Golub et al. (1999). Ribosomal protein S3 (RPS3, X57351) simultaneously participates in all three processes. This is reasonable since RPS3 is a housekeeping gene and ribosomal protein genes are usually over-expressed in some cancers. RPS3 encodes a ribosomal protein that is a component of the 40S subunit, where it forms a part of the domain where translation is initiated.

We used the Onto-Express (Draghici et al., 2003; Khatri et al., 2002) to investigate the enrichment of functional annotations of genes selected in each factor. Onto-Express starts by reading the input file that contains a list of GenBank accession numbers, and estimates the statistical significance of the enrichment of Gene-Ontology (GO) terms in the list with respect to a reference list. We used a list of all genes in the dataset as a reference array and hypergeometric distribution. Table 3 shows the enrichment of GO terms. We presented some significant biological processes for each factor, whose _P_-values were less than 0.01.

Table 3.

Enrichment of GO categories in genes selected by SNMF/R on the ALLAML leukemia dataset. We present some significant biological processes for each factor, whose _P_-values are less than 0.01

Table 3.

Enrichment of GO categories in genes selected by SNMF/R on the ALLAML leukemia dataset. We present some significant biological processes for each factor, whose _P_-values are less than 0.01

Figure 3 illustrates the matrices W and H obtained from SNMF/R with β = 1.0 on the CNS tumors dataset, which produced the lowest approximation error for five runs with different random initializations of W. Only one sample (the 18th sample, Brain_MGlio_8) was incorrectly assigned to the third cluster (rhabdoids). Our gene selection method suggested total of 367 genes (cluster1: 42, cluster2: 93, cluster3: 168, cluster4: 64). To more thoroughly characterize sets of genes dominantly expressed in different factors, we used the Onto-Express. The number of genes corresponding to each GO category was compared with the number of genes expected for the GO category in the Affymetrix HuGeneFL array. Significant differences from the expected were calculated with hypergeometric distribution. Table 4 shows biological processes with a significance of _P_-value < 0.01. The biological processes showing significant representations in the first factor were serotonin receptor signaling pathway, feeding behavior, and central nervous system development. Serotonin (5-hydroxytryptamine, or 5-HT) is a monoamine neurotransmitter and is known to regulate human mood, emotion, sleep and appetite in the central nervous system. Two GenBank accession numbers (U49516 and M81778) for serotonin receptor signaling pathway were linked to the same gene: 5-hydroxytryptamine (serotonin) receptor 2C (HTR2C). The GO category of feeding behavior seems to be related with childhood brain tumors known as medulloblastomas. Genes involved in the second factor were (4 genes: U52155, M81886, M64752, M81181) for potassium ion transport, (5 genes: X54673, M81886, M64752, M19650, L32961) for synaptic transmission, (3 genes: U62801, Z19002, M93426) for central nervous system development. This second cluster contains malignant glioma that is a tumor arising from glial cells. Genes corresponding to cell adhesion, cell motility and inflammatory response were highly expressed in the third factor. Genes highly expressed in normal cerebella were (5 genes: U79245, U33632, U90065, L36069, D79998) for potassium ion transport, (6 genes: M13577, U92457, L76627, U18244, M58583, U79667) for synaptic transmission, (3 genes: U52969, M13577, U76421) for central nervous system development and (6 genes: S81944, U79245, S95936, U33632, U90065, L36069) for ion transport. Detailed description of the clusters of samples and genes selected for each factor via SNMF/R can be found in supplementary materials. We have shown that SNMF/R can be used for clustering, cancer class discovery, gene selection and biological process analysis.

Table 4.

Enrichment of GO categories in genes selected by SNMF/R on the CNS tumors dataset. We present some significant biological processes for each factor, whose _P_-values are less than 0.01

Table 4.

Enrichment of GO categories in genes selected by SNMF/R on the CNS tumors dataset. We present some significant biological processes for each factor, whose _P_-values are less than 0.01

Fig. 3.

W (basis matrix) and H (coefficient matrix) obtained from SNMF/R with β = 1.0 for the CNS tumor dataset (34 samples: 10 classic medulloblastomas, 10 malignant gliomas, 10 rhabdoids, 4 normal cerebella) with the 5597 genes.

4 CONCLUSION

We present a novel sparse NMF algorithm via alternating non-negativity-constrained least squares. SNMF/R can be used for cancer class discovery and gene expression data analysis since it shows good biclustering performance and provides us with simple interpretation. This algorithm can be applied to many practical problems in bioinformatics and computational biology such as biomedical text mining and gene/protein microarray data analysis.

ACKNOWLEDGEMENTS

We would like to thank Dr Chris Ding, Dr Jean-Philippe Brunet, Dr Yuan Gao, Prof. Lars Eldén and Prof. Robert J. Plemmons for their valuable comments. In particular, we would also like to thank Prof. Chih-Jen Lin, Prof. Paul Tseng and Prof. Luigi Grippo for discussions on the convergence property. This material is based upon work supported in part by the National Science Foundation Grants ACI-0305543 and CCF-0621889. Any opinions, findings and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the National Science Foundation.

Conflict of Interest: none declared.

REFERENCES

et al.

Matrices, vector spaces, and information retrieval

,

SIAM Rev

,

1999

, vol.

41

(pg.

335

-

362

)

et al.

Algorithms and applications for approximate nonnegative matrix factorization

,

Comput. Stat. Data Anal

,

2006

(to appear)

A fast non-negativity-constrained least squares algorithm

,

J. Chemometrics

,

1997

, vol.

11

(pg.

393

-

401

)

et al.

Metagenes and molecular pattern discovery using matrix factorization

,

Proc. Natl Acad. Sci. USA

,

2004

, vol.

101

(pg.

4164

-

4169

)

et al.

Biclustering of gene expression data by non-smooth non-negative matrix factorization

,

BMC Bioinformatics

,

2006

, vol.

7

pg.

78

et al.

High-resolution genomic profiles define distinct clinico-pathogenetic subgroups of multiple myeloma patients

,

Cancer Cell

,

2006

, vol.

9

(pg.

313

-

325

)

et al.

Discovering semantic features in the literature: a foundation for building functional associations

,

BMC Bioinformatics

,

2006

, vol.

7

pg.

41

et al.

Adaptive dimension reduction for clustering high dimensional data

,

2002

In Proceedings of the 2nd IEEE International Conference on Data Mining

Maebashi, Japan

et al.

Onto-tools, the toolkit of the modern biologist: Onto-Express, Onto-Compare, Onto-Design and Onto-Translate

,

Nucleic Acids Res

,

2003

, vol.

31

(pg.

3775

-

3781

)

et al.

Multi-way clustering of microarray data using probabilistic sparse matrix factorization

,

Bioinformatics

,

2005

, vol.

21

Suppl. 1

(pg.

i144

-

i151

)

Improving molecular cancer class discovery through sparse non-negative matrix factorization

,

Bioinformatics

,

2005

, vol.

21

(pg.

3970

-

3975

)

et al.

Molecular classification of cancer: class discovery and class prediction by gene expression monitoring

,

Science

,

1999

, vol.

286

(pg.

531

-

537

)

Accelerating the Lee-Seung algorithm for non-negative matrix factorization

,

Technical report

,

2005

Department of Computational and Applied Mathematics, Rice University

On the convergence of the block nonlinear Gauss-Seidel method under convex constraints

,

Operations Res. Lett

,

2000

, vol.

26

(pg.

127

-

136

)

Non-negative matrix factorization with sparseness constraints

,

J. Machine Learning Res

,

2004

, vol.

5

(pg.

1457

-

1469

)

et al.

Profiling gene expression using onto-express

,

Genomics

,

2002

, vol.

79

(pg.

266

-

270

)

et al.

Dimension reduction in text classification with support vector machines

,

J. Machine Learning Res

,

2005

, vol.

6

(pg.

37

-

53

)

Subsystem identification through dimensionality reduction of large-scale gene expression data

,

Genome Res

,

2003

, vol.

13

(pg.

1706

-

1718

)

,

Solving Least Squares Problems

,

1974

NJ

Prentice-Hall, Englewood Cliffs

Learning the parts of objects by non-negative matrix factorization

,

Nature

,

1999

, vol.

401

(pg.

788

-

791

)

Algorithms for non-negative matrix factorization

,

2000

Proceedings of Neural Information Processing Systems

(pg.

556

-

562

)

et al.

Learning spatially localized parts-based representations

,

2001

In Proceedings IEEE Conference on Computer Vision and Pattern Recognition (CVPR)

(pg.

207

-

212

)

et al.

Marked genomic differences characterize primary and secondary glioblastoma subtypes and identify two distinct molecular and clinical secondary glioblastoma entities

,

Cancer Res

,

2006

, vol.

66

(pg.

11502

-

11513

)

MATLAB

,

User's Guide.

,

1992

The MathWorks, Inc

Natick, MA 01760

et al.

Nonsmooth nonnegative matrix factorization (nsNMF)

,

IEEE, Trans. Pattern Anal. Machine Intell

,

2006

, vol.

28

(pg.

403

-

415

)

et al.

Text mining using non-negative matrix factorizations

,

2004

Proceedings SIAM International Conference on Data Mining (SDM'04)

et al.

Nonnegative matrix factorization for spectral data analysis

,

Linear Algebra and Applications

,

2006

(to appear)

et al.

Theme discovery from gene lists for identification and viewing of multiple functional groups

,

BMC Bioinformatics

,

2005

, vol.

6

pg.

162

et al.

Prediction of central nervous system embryonal tumour outcome based on gene expression

,

Nature

,

2002

, vol.

415

(pg.

436

-

442

)

Regression shrinkage and selection via LASSO

,

J. Roy. Statist. Soc. B

,

1996

, vol.

58

(pg.

267

-

288

)

Fast algorithm for the solution of large-scale non-negativity-constrained least squares problems

,

J. Chemometrics

,

2004

, vol.

18

(pg.

441

-

450

)

Author notes

Associate Editor: David Rocke

© The Author 2007. Published by Oxford University Press. All rights reserved. For Permissions, please email: journals.permissions@oxfordjournals.org

Citations

Views

Altmetric

Metrics

Total Views 22,052

18,322 Pageviews

3,730 PDF Downloads

Since 11/1/2016

Citations

599 Web of Science

×

Email alerts

Citing articles via

• Latest

• Most Read

• Most Cited

More from Oxford Academic