Non-negative matrix factorization with sparseness constraints

Patrik O. Hoyer
HIIT Basic Research Unit
Department of Computer Science
University of Helsinki, Finland
www.cs.helsinki.fi/patrik.hoyer

Abstract

Non-negative matrix factorization (NMF) is a recently developed technique for finding parts-based, linear representations of non-negative data. Although it has successfully been applied in several applications, it does not always result in parts-based representations. In this paper, we show how explicitly incorporating the notion of 'sparseness' improves the found decompositions. Additionally, we provide complete MATLAB code both for standard NMF and for our extension. Our hope is that this will further the application of these methods to solving novel data-analysis problems.

Keywords: Non-negative matrix factorization, sparseness, data-adaptive representations

Executive Summary: Non-negative matrix factorization (NMF) extracts additive, parts-based representations from non-negative data matrices but often fails to deliver intuitive decompositions because it provides no direct control over how sparse those representations become. This paper develops an extension that adds explicit sparseness constraints on the basis vectors or coefficients (or both), together with a complete MATLAB implementation of the method and its predecessors.

The work evaluates whether enforcing chosen sparseness levels improves the quality of the learned factors on image data where standard NMF produces global or otherwise unhelpful features. The approach uses a projected gradient algorithm that alternates multiplicative updates with a projection step enforcing non-negativity plus prescribed L1 and L2 norms; the projection itself is shown to converge in a handful of iterations even at dimensionalities up to 10,000. Experiments were run on the CBCL and ORL face databases and on high-pass-filtered natural-image patches split into ON and OFF contrast channels.

On the ORL faces, raising the sparseness of the basis vectors from roughly 0.5 to 0.75 switched the representation from global to clearly local (eyes, noses, mouths). On the natural-image data, constraining only the coefficients to a sparseness of 0.85 produced oriented, Gabor-like filters instead of the circular blobs returned by ordinary NMF. The same framework can also be tuned to recover global representations when desired. Relative to earlier extensions such as LNMF and SNMF, the new method supplies explicit rather than implicit sparseness control and succeeds on a broader range of data sets.

These results indicate that analysts can now obtain more interpretable, task-appropriate factors without extensive trial-and-error parameter search. The accompanying software removes a practical barrier to adoption by supplying ready-to-run code for all compared variants.

The principal limitations are that suitable sparseness targets must still be chosen by the user and that the reported gains rest on three image collections; performance on other data types remains to be verified. Overall is high for the image-domain findings and for the algorithmic claims, moderate for broader applicability.

1. Introduction

Section Summary: The introduction explains that many data analysis tasks require finding compact representations that reveal hidden structure in the data. Non-negative matrix factorization achieves this by breaking a data matrix into additive, non-negative parts, which often produces useful sparse encodings, but without any way to adjust how sparse those encodings are. The paper therefore extends the method to let users control sparseness directly, shows that this yields clearer parts-based representations, and supplies accompanying MATLAB software.

A fundamental problem in many data-analysis tasks is to find a suitable representation of the data. A useful representation typically makes latent structure in the data explicit, and often reduces the dimensionality of the data so that further computational methods can be applied.

Non-negative matrix factorization (NMF) [1, 2] is a recent method for finding such a representation. Given a non-negative data matrix ${V}$, NMF finds an approximate factorization ${V} \approx {W} {H}$ into non-negative factors ${W}$ and ${H}$. The non-negativity constraints make the representation purely additive (allowing no subtractions), in contrast to many other linear representations such as principal component analysis (PCA) and independent component analysis (ICA) [3].

One of the most useful properties of NMF is that it usually produces a sparse representation of the data. Such a representation encodes much of the data using few 'active' components, which makes the encoding easy to interpret. Sparse coding [4] has also, on theoretical grounds, been shown to be a useful middle ground between completely distributed representations, on the one hand, and unary representations (grandmother cells) on the other [5, 6]. However, because the sparseness given by NMF is somewhat of a side-effect rather than a goal, one cannot in any way control the degree to which the representation is sparse. In many applications, more direct control over the properties of the representation is needed.

In this paper, we extend NMF to include the option to control sparseness explicitly. We show that this allows us to discover parts-based representations that are qualitatively better than those given by basic NMF. We also discuss the relationship between our method and other recent extensions of NMF [7, 8, 9].

Additionally, this contribution includes a complete MATLAB package for performing NMF and its various extensions. Although the most basic version of NMF requires only two lines of code and certainly does not warrant distributing a separate software package, its several extensions involve more complicated operations; the absense of ready-made code has probably hindered their widespread use so far. We hope that our software package will alleviate the problem.

This paper is structured as follows. In Section 2 we describe non-negative matrix factorization, and discuss its success but also its limitations. Section 3 discusses why and how to incorporate sparseness constraints into the NMF formulation. Section 4 provides experimental results that verify our approach. Finally, Section 5 and Section 6 compare our approach to other recent extensions of NMF and conclude the paper.

2. Non-negative matrix factorization

Section Summary: Non-negative matrix factorization approximates collections of non-negative data measurements as additive combinations of a smaller set of non-negative basis vectors, expressed by factoring the data matrix V into the product of two non-negative matrices W and H that minimize reconstruction error. Unlike principal component analysis or independent component analysis, the non-negativity constraint means the data can only be built from parts through addition, never subtraction, which often produces intuitive parts-based features such as facial elements in image datasets. The approach succeeds on some real-world data like faces or documents but can yield less intuitive global or blob-like results on others, depending on factors such as alignment or data preprocessing.

Non-negative matrix factorization is a linear, non-negative approximate data representation. Let's assume that our data consists of $T$ measurements of $N$ non-negative scalar variables. Denoting the ($N$-dimensional) measurement vectors ${v}^t$ ($t=1,\ldots, T$), a linear approximation of the data is given by

$ {v}^t \approx \sum_{i=1}^M {w}_i h_i^t = {W} {h}^t, $

where ${W}$ is an $N\times M$ matrix containing the basis vectors ${w}_i$ as its columns. Note that each measurement vector is written in terms of the same basis vectors. The $M$ basis vectors ${w}_i$ can be thought of as the 'building blocks' of the data, and the ($M$-dimensional) coefficient vector ${h}^t$ describes how strongly each building block is present in the measurement vector ${v}^t$.

Arranging the measurement vectors ${v}^t$ into the columns of an $N\times T$ matrix ${V}$, we can now write

$ {V} \approx {W} {H}, $

where each column of ${H}$ contains the coefficient vector ${h}^t$ corresponding to the measurement vector ${v}^t$. Written in this form, it becomes apparent that a linear data representation is simpy a factorization of the data matrix. Principal component analysis, independent component analysis, vector quantization, and non-negative matrix factorization can all be seen as matrix factorization, with different choices of objective function and/or constraints.

Whereas PCA and ICA do not in any way restrict the signs of the entries of ${W}$ and ${H}$, NMF requires all entries of both matrices to be non-negative. What this means is that the data is described by using additive components only. This constraint has been motivated in a couple of ways: First, in many applications one knows (e.g. by the rules of physics) that the quantities involved cannot be negative. In such cases, it can be difficult to interpret the results of PCA and ICA [1, 10]. Second, non-negativity has been argued for based on the intuition that parts are generally combined additively (and not subtracted) to form a whole; hence, these constraints might be useful for learning parts-based representations [2].

Given a data matrix ${V}$, the optimal choice of matrices ${W}$ and ${H}$ are defined to be those non-negative matrices that minimize the reconstruction error between ${V}$ and ${W} {H}$. Various error functions have been proposed [1, 11], perhaps the most widely used is the squared error (euclidean distance) function

$ E({W}, {H}) = | {V} - {W} {H}|^2 = \sum_{i,j} (V_{ij} - ({W} {H})_{ij})^2. $

Although the minimization problem is convex in ${W}$ and ${H}$ separately, it is not convex in both simultaneously. Paatero [1] gave a gradient algorithm for this optimization, whereas Lee and Seung [11] devised a multiplicative algorithm that is somewhat simpler to implement and also showed good performance.

Although some theoretical work on the properties of the NMF representation exists [15], much of the appeal of NMF comes from its empirical success in learning meaningful features from a diverse collection of real-life datasets. Lee and Seung [2] showed that, when the dataset consisted of a collection of face images [12], the representation consisted of basis vectors encoding for the mouth, nose, eyes, etc; the intuitive features of face images. In Figure 1 a we have reproduced that basic result using the same dataset. Additionally, they showed that meaningful topics can be learned when text documents are used as data. Subsequently, NMF has been successfully applied to a variety of datasets [16, 17, 18, 19].

Despite this success, there also exist datasets for which NMF does not give an intuitive decomposition into parts that would correspond to our idea of the 'building blocks' of the data. In [7], the authors showed that when NMF was applied to a different facial image database [13], the representation was global rather than local, qualitatively different from that reported in [2]. Again, we have rerun that experiment and confirm those results, see Figure 1 b. The difference was mainly attributed to how well the images were hand-aligned [7].

Another case where the decomposition found by NMF does not match the underlying elements of the data is shown in Figure 1 c. In this experiment [14], natural image patches were high-pass filtered and subsequently split into positive ('ON') and negative ('OFF') contrast channels, in a process similar to how visual information is processed by the retina. When NMF is applied to such a dataset, the resulting decomposition does not consist of the oriented filters which form the cornerstone of most of modern image processing. Rather, NMF represents these images using simple, dull, circular 'blobs'.

We will show that, in both of the above cases, explicitly controlling the sparseness of the representation leads to representations that are parts-based and match the intuitive features of the data.

3. Adding sparseness constraints to NMF

Section Summary: The section explains sparseness as a property in which most elements of a vector stay near zero while only a few take on large values, and it introduces a simple numerical measure based on the relationship between the L1 and L2 norms. It then extends non-negative matrix factorization so that the user can optionally force the basis matrix, the coefficient matrix, or both to meet chosen sparseness levels while still minimizing reconstruction error. The resulting method uses a projected gradient algorithm that repeatedly adjusts the factors and projects them onto the required sparseness constraints.

3.1 Sparseness

The concept of 'sparse coding' refers to a representational scheme where only a few units (out of a large population) are effectively used to represent typical data vectors [4]. In effect, this implies most units taking values close to zero while only few take significantly non-zero values. Figure 2 illustrates the concept and our sparseness measure (defined below).

Numerous sparseness measures have been proposed and used in the literature to date. Such measures are mappings from $\mathbb{R}^n$ to $\mathbb{R}$ which quantify how much energy of a vector is packed into only a few components. On a normalized scale, the sparsest possible vector (only a single component is non-zero) should have a sparseness of one, whereas a vector with all elements equal should have a sparseness of zero.

In this paper, we use a sparseness measure based on the relationship between the $L_1$ norm and the $L_2$ norm:

$ \text{sparseness}({x}) = \frac{\sqrt{n} - \left(\sum |x_i|\right)/\sqrt{\sum x_i^2}}{ \sqrt{n} - 1 }, $

where $n$ is the dimensionality of ${x}$. This function evaluates to unity if and only if ${x}$ contains only a single non-zero component, and takes a value of zero if and only if all components are equal (up to signs), interpolating smoothly between the two extremes.

3.2 NMF with sparseness constraints

Our aim is to constrain NMF to find solutions with desired degrees of sparseness. The first question to answer is then: what exactly should be sparse? The basis vectors ${W}$ or the coefficients ${H}$? This is a question that cannot be given a general answer; it all depends on the specific application in question. Further, just transposing the data matrix switches the role of the two, so it is easy to see that the choice of which to constrain (or both, or none) must be made by the experimenter.

For example, a doctor analyzing disease patterns might assume that most diseases are rare (hence sparse) but that each disease can cause a large number of symptoms. Assuming that symptoms make up the rows of her matrix and the columns denote different individuals, in this case it is the 'coefficients' which should be sparse and the 'basis vectors' unconstrained. On the other hand, when trying to learn useful features from a database of images, it might make sense to require both ${W}$ and ${H}$ to be sparse, signifying that any given object is present in few images and affects only a small part of the image.

These considerations lead us to defining NMF with sparseness constraints as follows:

{{Definition:} NMF with sparseness constraints

Given a non-negative data matrix ${V}$ of size $N \times T$, find the non-negative matrices ${W}$ and ${H}$ of sizes $N \times M$ and $M \times T$ (respectively) such that

$ E({W}, {H}) = | {V} - {W} {H}|^2 $

is minimized, under optional constraints

$ \begin{aligned} \text{sparseness}({w}_i) & = S_w, ;; \forall i \ \text{sparseness}({h}_i) & = S_h, ;; \forall i, \end{aligned} $

where ${w}_i$ is the $i$:th column of ${W}$ and ${h}_i$ is the $i$:th row of ${H}$. Here, $M$ denotes the number of components, and $S_w$ and $S_h$ are the desired sparsenesses of ${W}$ and ${H}$ (respectively). These three parameters are set by the user.

Note that we did not constrain the scales of ${w}_i$ or ${h}_i$ yet. However, since ${w}_i {h}_i = ({w}_i\lambda)({h}_i/\lambda)$ we are free to arbitrarily fix any norm of either one. In our algorithm, we thus choose to fix the $L_2$ norm of ${h}_i$ to unity, as a matter of convenience.

3.3 Algorithm

We have devised a projected gradient descent algorithm for NMF with sparseness constraints. This algorithm essentially takes a step in the direction of the negative gradient, and subsequently projects onto the constraint space, making sure that the taken step is small enough that the objective function Equation (5) is reduced at every step. The main muscle of the algorithm is the projection operator which enforces the desired degree of sparseness. This operator is described in detail following this algorithm.

{{Algorithm:} NMF with sparseness constraints

::: {.visual-block}

:::

Above, $\otimes$ and $\oslash$ denote elementwise multiplication and division, respectively. Moreover, $\mu_{{W}}$ and $\mu_{{H}}$ are small positive constants (stepsizes) which must be set appropriately for the algorithm to work. Fortunately, they need not be set by the user; our implementation of the algorithm automatically adapts these parameters. The multiplicative steps are directly taken from [11] and are used when constraints are not to be applied.

Many of the steps in the above algorithm require a projection operator which enforces sparseness by explicitly setting both $L_1$ and $L_2$ norms (and enforcing non-negativity). This operator is defined as follows

::: {.visual-block}

:::

In words, the above algorithm works as follows: We start by projecting the given vector onto the hyperplane $\sum s_i = L_1$. Next, within this space, we project to the closest point on the joint constraint hypersphere (intersection of the sum and the $L_2$ constraints). This is done by moving radially outward from the center of the sphere (the center is given by the point where all components have equal values). If the result is completely non-negative, we have arrived at our destination. If not, those components that attained negative values must be fixed at zero, and a new point found in a similar fashion under those additional constraints.

Note that, once we have a solution to the above non-negative problem, it would be straightforward to extend it to a general solution without non-negativity constraints. If a given component of ${x}$ is positive (negative), we know because of the symmetries of $L_1$ and $L_2$ norms that the optimal solution ${s}$ will have the corresponding component positive or zero (negative or zero). Thus, we may simply record the signs of ${x}$, take the absolute value, perform the projection in the first quadrant using the algorithm above, and re-enter the signs into the solution.

In principle, the devised projection algorithm may take as many as $\text{dim}({x})$ iterations to converge to the correct solution (because at each iteration the algorithm either converges, or at least one component is added to the set of zero valued components). In practice, however, the algorithm converges much faster. In Section 4 we show that even for extremely high dimensions the algorithm typically converges in only a few iterations.

3.4 Matlab implementation

Our software package, available at http://www.cs.helsinki.fi/patrik.hoyer/ implements all the details of the above algorithm. In particular, we monitor the objective function $E$ throughout the optimization, and adapt the stepsizes to ensure convergence. The software package contains, in addition to the projection operator and NMF code, all the files needed to reproduce the results described in this paper, with the exception of datasets. For copyright reasons the face image databases are not included, but they can easily be downloaded separately from their respective www addresses.

4. Experiments with sparseness constraints

Section Summary: Adding sparseness constraints to NMF enables it to discover localized, parts-based features from data where standard NMF instead produces global or circular patterns, as shown on misaligned face images and natural image patches. By adjusting the required sparseness level separately on the basis vectors or the coefficients, the method can flexibly produce either local edge-like features or more holistic representations. Experiments also confirm that the projection step used to enforce these constraints converges in fewer than ten iterations even in high-dimensional cases.

In this section, we show that adding sparseness constraints to NMF can make it find parts-based representations in cases where unconstrained NMF does not. In addition, we experimentally verify our claim that the projection operator described in Section 3.3 converges in only a few iterations even when the dimensionality of the vector is high.

4.1 Representations learned from face image databases

Recall from Section 2 the mixed results of applying standard NMF to face image data. Lee and Seung [2] originally showed that NMF found a parts-based representation when trained on data from the CBCL database. However, when applied to the ORL dataset, in which images are not as well aligned, a global decomposition emerges. These results were shown in Figure 1 a and Figure 1 b. To compare, we applied sparseness constrained NMF to both face image datasets.

For the CBCL data, some resulting bases are shown in Figure 3. Setting a high sparseness value for the basis images results in a local representation similar to that found by standard NMF. However, we want to emphasize the fact that sparseness constrained NMF does not always lead to local solutions: Global solutions can be obtained by deliberately setting a low sparseness on the basis images, or by requiring a high sparseness on the coefficients (forcing each coefficient to try to represent more of the image).

The ORL database provides the more interesting test of the method. In Figure 4 we show bases learned by sparseness constrained NMF, for various sparseness settings. Note that our method can learn a parts-based representation of this dataset, in contrast to standard NMF. Also note that the representation is not very sensitive to the specific sparseness level chosen.

4.2 Basis derived from natural image patches

In Figure 1 c we showed that standard NMF applied to natural image data produces only circular features, not oriented features like those employed by modern image processing techniques. Here, we tested the result of using additional sparseness constraints. Figure 5 shows the basis vectors obtained by putting a sparseness constraint on the coefficients ($S_h = 0.85$) but leaving the sparseness of the basis vectors unconstrained. In this case, NMF learns oriented features that represent edges and lines. Such oriented features are widely regarded as the best type of low-level features for representing natural images, and similar features are also used by the early visual system of the biological brain [20, 21, 22, 23]. This example illustrates that sparseness constrained NMF does not simply 'sparsify' the result of standard, unconstrained NMF, but rather can find qualitatively different parts-based representations that are more compatible with the sparseness assumptions.

4.3 Convergence of algorithm implementing the projection step

To verify the performance of our projection method we performed extensive tests, varying the number of dimensions, the desired degree of sparseness, and the sparseness of the original vector. The desired and the initial degrees of sparseness were set to 0.1, 0.3, 0.5, 0.7, and 0.9, and the dimensionality of the problem was set to 2, 3, 5, 10, 50, 100, 500, 1000, 3000, 5000, and 10000. All combinations of sparsenesses and dimensionalities were analyzed. Based on this analysis, the worst case (most iterations on average required) was when the desired degree of sparseness was high (0.9) but the initial sparseness was low (0.1). In Figure 6 we plots the number of iterations required for this worst case, as a function of dimensionality. Even in this worst-case scenario, and even for the highest tested dimensionality, the algorithm never required more than 10 iterations to converge. Thus, although we do not have analytical bounds on the performance on the algorithm, empirically the projection method performs extremely well.

5. Relation to other recent work

Section Summary: This section compares the current non-negative matrix factorization approach to prior extensions of NMF, such as local NMF and sparse NMF, noting that those methods achieve local or sparse features only indirectly and often fail to produce oriented components from natural images. It highlights the present model's advantage in providing explicit, adjustable control over sparseness without trial-and-error tuning, along with shared code for direct comparisons. The work also relates closely to independent component analysis and its non-negative variants, differing mainly in how signs are restricted and how sparseness is enforced.

5.1 Extensions of NMF

Several authors have noted the shortcomings of standard NMF, and suggested extensions and modifications of the original model. Li et al [7] noted that NMF found only global features from the ORL database (see Figure 1 b) and suggested an extension they call Local Non-negative Matrix Factorization (LNMF). Their method indeed produces local features from the ORL database, similar to those given by our method (Figure 4 c). However, it does not produce oriented filters from natural image data (results not shown). Further, there is no way to explicitly control the sparseness of the representation, should this be needed.

In [8] the current author extended the NMF framework to include an adjustable sparseness parameter. The present paper is an extension of those ideas. The main improvement is that in the present model sparseness is adjusted explicitly, rather than implicitly. This means that one does not any more need to employ trial-and-error to find the parameter setting that yields the desired level of sparseness.

Finally, Liu et al [9] also noted the need for incorporating the notion of sparseness, and suggested an extension termed Sparse Non-negative Matrix Factorization (SNMF). Their extension is similar in spirit and form to that given in [8] with the added benefit of yielding a more convenient, faster algorithm. Nevertheless, it also suffers from the drawback that sparseness is only controlled implicitly. Furthermore, their method does not yield oriented features from natural image data (results not shown).

In summary, the framework presented in the present paper improves on these previous extensions by allowing explicit control of the statistical properties of the representation.

In order to facilitate the use of, and comparison between, the various extensions of NMF, they are all provided as part of the Matlab code package distributed with this paper. Using this package readers can effortlessly verify our current claims by applying the algorithms to the various datasets. Moreover, the methods can be compared head-to-head on new interesting datasets.

5.2 Non-negative independent component analysis

Our method has a close connection to the statistical technique called independent component analysis (ICA) [3]. ICA attempts to find a matrix factorization similar to ours, but with two important differences. First, the signs of the components are in general not restricted; in fact, symmetry is often assumed, implying an approximately equal number of positive and negative elements. Second, the sources are not forced to any desired degree of sparseness (as in our method) but rather sparseness is incorporated into the objective function to be optimized. The sparseness goal can be put on either ${W}$ or ${H}$, or both [24].

Recently, some authors have considered estimating the ICA model in the case of one-sided, non-negative sources [25, 26]. In these methods, non-negativity is not specified as a constraint but rather as an objective; hence, complete non-negativity of the representation is seldom achieved for real-life datasets. Nevertheless, one can show that if the linear ICA model holds, with non-negative components, these methods can identify the model.

6. Conclusions

Section Summary: Non-negative matrix factorization, or NMF, is a helpful method for breaking down various types of data into understandable parts because it often produces sparse results that are easy to interpret. However, when standard NMF does not provide enough sparsity, researchers may need more direct control over how sparse the results are. This paper addresses that need by introducing a mathematical tool to enforce exact levels of sparsity in NMF models, demonstrating its ability to create new kinds of data representations, and releasing software to make these techniques accessible to others in the field.

Non-negative matrix factorization (NMF) has proven itself a useful tool in the analysis of a diverse range of data. One of its most useful properties is that the resulting decompositions are often intuitive and easy to interpret because they are sparse. Sometimes, however, the sparseness achieved by NMF is not enough; in such situations it might be useful to control the degree of sparseness explicitly. Our main contributions of this paper were (a) to describe a projection operator capable of simultaneously enforcing both $L_1$ and $L_2$ norms and hence any desired degree of sparseness, (b) to show its use in the NMF framework for learning representations that could not be obtained by regular NMF, and (c) to provide a software package to enable researchers and practitioners to easily perform NMF and its various extensions. We hope that all three contributions will prove useful to the field of data-analysis.

Acknowledgements

The author wishes to thank Jarmo Hurri, Aapo Hyvärinen, and Fabian Theis for useful discussions and comments on the manuscript.

References

Section Summary: The references section compiles 26 academic sources, primarily journal articles, conference papers, and books from the 1980s through 2004, that underpin research on non-negative matrix factorization and related decomposition methods. These works address core techniques such as sparse coding, independent component analysis, and parts-based representations, along with their uses in areas like image processing, neuroscience, and gene data analysis. They feature foundational contributions from researchers including Lee and Seung, Hyvärinen, and Olshausen, plus links to public datasets for faces and other visual stimuli.

[1] P. Paatero and U. Tapper. Positive matrix factorization: A non-negative factor model with optimal utilization of error estimates of data values. Environmetrics, 5:111–126, 1994.

[2] D. D. Lee and H. S. Seung. Learning the parts of objects by non-negative matrix factorization. Nature, 401(6755):788–791, 1999.

[3] A. Hyvärinen, J. Karhunen, and E. Oja. Independent Component Analysis. Wiley Interscience, 2001.

[4] D. J. Field. What is the goal of sensory coding? Neural Computation, 6:559–601, 1994.

[5] P. Földiák and M. P. Young. Sparse coding in the primate cortex. In Michael A. Arbib, editor, The Handbook of Brain Theory and Neural Networks, pages 895–898. The MIT Press, Cambridge, Massachusetts, 1995.

[6] S. Thorpe. Localized versus distributed representations. In Michael A. Arbib, editor, The Handbook of Brain Theory and Neural Networks, pages 549–552. The MIT Press, Cambridge, Massachusetts, 1995.

[7] S. Z. Li, X. Hou, H. Zhang, and Q. Cheng. Learning spatially localized parts-based representations. In Proc. IEEE Conf. on Computer Vision and Pattern Recognition (CVPR), Vol. I, pages 207–212, Hawaii, USA, 2001.

[8] P. O. Hoyer. Non-negative sparse coding. In Neural Networks for Signal Processing XII (Proc. IEEE Workshop on Neural Networks for Signal Processing), pages 557–565, Martigny, Switzerland, 2002.

[9] W. Liu, N. Zheng, and X. Lu. Non-negative matrix factorization for visual coding. In Proc. IEEE Int. Conf. on Acoustics, Speech and Signal Processing (ICASSP'2003), 2003.

[10] L. Parra, C. Spence, P. Sajda, A. Ziehe, and K.-R. Müller. Unmixing hyperspectral data. In Advances in Neural Information Processing 12 (Proc. NIPS*99), pages 942–948. MIT Press, 2000.

[11] D. D. Lee and H. S. Seung. Algorithms for non-negative matrix factorization. In Advances in Neural Information Processing 13 (Proc. NIPS*2000). MIT Press, 2001.

[12] MIT Center For Biological and Computation Learning. CBCL face database #1. Available at: http://cbcl.mit.edu/cbcl/software-datasets/FaceData2.html.

[13] AT&T Laboratories Cambridge. The ORL database of faces. Available at: http://www.uk.research.att.com/facedatabase.html.

[14] P. O. Hoyer. Modeling receptive fields with non-negative sparse coding. In E. De Schutter, editor, Computational Neuroscience: Trends in Research 2003. Elsevier, Amsterdam, 2003. Also published in: Neurocomputing 52-54 (2003), pp 547-552.

[15] D. Donoho and V. Stodden. When does non-negative matrix factorization give a correct decomposition into parts? In Advances in Neural Information Processing 16 (Proc. NIPS*2003). MIT Press, 2004.

[16] G. Buchsbaum and O. Bloch. Color categories revealed by non-negative matrix factorization of munsell color spectra. Vision Research, 42:559–563, 2002.

[17] J. P. Brunet, P. Tamayo, T. R. Golub, and J. P. Mesirov. Metagenes and molecular pattern discovery using matrix factorization. Proceedings of the National Academy of Sciences, 101:4164–4169, 2004.

[18] K. Jung and E. Y. Kim. Automatic text extraction for content-based image indexing. Advances in Knowledge Discovery and Data Mining: Proceedings Lecture Notes in Artificial Intelligence, 3056:497–507, 2004.

[19] P. M. Kim and B. Tidor. Subsystem identification through dimensionality reduction of large-scale gene expression data. Genome Research, 13:1706–1718, 2003.

[20] D. J. Field. Relations between the statistics of natural images and the response properties of cortical cells. Journal of the Optical Society of America, 4:2379–2394, 1987.

[21] E. P. Simoncelli, W. T. Freeman, E. H. Adelson, and D. J. Heeger. Shiftable multiscale transforms. IEEE Transactions on Information Theory, 38:587–607, 1992.

[22] B. A. Olshausen and D. J. Field. Emergence of simple-cell receptive field properties by learning a sparse code for natural images. Nature, 381:607–609, 1996.

[23] A. J. Bell and T. J. Sejnowski. The 'independent components' of natural scenes are edge filters. Vision Research, 37:3327–3338, 1997.

[24] J. V. Stone, J. Porrill, N. R. Porter, and I. D. Wilkinson. Spatiotemporal independent component analysis of event-related fMRI data using skewed probability density functions. Neuroimage, 15:407–421, 2002.

[25] M. Plumbley. Algorithms for non-negative independent component analysis. IEEE Transactions on Neural Networks, 14(3):534–543, 2003.

[26] E. Oja and M. Plumbley. Blind separation of positive sources by globally convergent gradient search. Neural Computation, 16(9):1811–1825, 2004.