Fast approximations to structured sparse coding and applications to object classification

Arthur Szlam, Karol Gregor, Yann LeCun

Abstract

We describe a method for fast approximation of sparse coding. The input space is subdivided by a binary decision tree, and we simultaneously learn a dictionary and assignment of allowed dictionary elements for each leaf of the tree. We store a lookup table with the assignments and the pseudoinverses for each node, allowing for very fast inference. %During inference %the algorithm uses a tree to assign an input to a group of allowed dictionary elements and then finds the corresponding coefficient values %using a cached pseudoinverse. We give an algorithm for learning the tree, the dictionary and the dictionary element assignment, and In the process of describing this algorithm, we discuss the more general problem of learning the groups in group structured sparse modelling. We show that our method creates good sparse representations by using it in the object recognition framework of [lazebnik06,yang-cvpr-09]. Implementing our own fast version of the SIFT descriptor the whole system runs at $20$ frames per second on $321 \times 481$ sized images on a laptop with a quad-core cpu, while sacrificing very little accuracy on the Caltech 101 and 15 scenes benchmarks.

Introduction

Sparse modeling [3, 4] has proven to be a useful framework for signal processing. Each point from a dataset consisting of vectors in a Euclidean space is represented by a vector with only a few nonzero coefficients. Sparse modeling has lead to state of the art algorithms in image denoising, inpainting, supervised learning, and of particular interest here, object recognition. The systems described in [1, 2, 5, 6, 7] use sparse coding as an integral element. Since the coding is done densely in an image with relatively large dictionaries, this is a computationally expensive part of the recognition system, and a barrier to real time application. The main contribution of this paper is a

(i.e. sets of nonzero coordinates of z ). This can be a problem in classification tasks. Other times we may know in advance some structure in the data that the coefficients should preserve. Various forms of structured sparsity are explored in [8, 9, 10, 11].

A simple form of structured sparsity is given by specifying a list of L allowable active sets, and some function g : R d ↦→ 1 , ..., L associating to each x to one of the L configurations. An example of this is the output of many subspace clustering algorithms. There, X is reordered and partitioned into PX = [ X 1 X 2 ... X L ] (where P is a permutation matrix), so that each block X j is near a low dimensional subspace spanned by B j . Supposing for simplicity that each of the B j are of the same dimension q , then if we set W = [ B 1 ...B L ] , the allowable active sets are given by 1 , ..., q , q + 1 , ..., 2 q , etc. By setting the allowable active sets to the blocks, and the function g to simply map each point to its nearest subspace (say in the standard sense of Euclidean projections), then we get an example of structured sparsity as described above; this sort of method is used in object recognition in [6].

In this work we will try to learn the L configurations as well as the dictionary. We introduce a LLoydlike algorithm that alternates between updating the dictionary, updating the assignments of each data point to the groups, and updating the dictionary elements associated to a group via simultaneous OMP [12].

At inference time, we need a fast method for determining which group an x belongs to. This is computationally expensive if there is a large number of groups and one needs check the projection onto each group. However, by specializing the Lloyd type algorithm to the case when each group is composed of a union of (perhaps only one) leaves of a binary decision tree, we will build a fast inference scheme into the learned dictionary. The key idea is that by using SOMP, we can learn which leaves should use which dictionary elements as we train the dictionary. To code an input, we march it down the tree until we arrive at the appropriate leaf. In addition to the decision vectors and thresholds, we will store a lookup table with the active set of each leaf as learned above, and the pseudoinverse of the columns of W corresponding to that active set. Thus after following x down the tree we need only make one matrix multiplication to get the coefficients.

Finally, we would like use these algorithms to build an accurate real time recognition system. We focus on a par- ticular architecture studied in [1, 2, 6, 7]. First, SIFT descriptors are calculated densely over the image. Then (a form of) sparse coding is used to calculate a sparse vector at every location from the corresponding sift vector. Then each feature is pooled over a small number of spatial regions and the results are concatenated. Finally the labels are obtained using linear SVM or logistic regression.

We use this pipeline with two modifications. First we write our own fast implementation of the SIFT descriptor. Second we use our fast algorithm for the sparse coding step. The resulting system achieves nearly the same performance as exact sparse coding calculation but processes 321 × 481 size images at the rate of 20 frames per second on a laptop computer with a quad core cpu.

The rest of this paper is organized as follows: Section 2, we discuss greedy structured sparse modeling, and describe in depth how to train a model that learns the structure, and that respects a given set of groups given by a tree. In section 3, we show experiments on image patches to qualitatively demonstrate what learned groups look like and then we apply our methods to object recognition.

Here we will first suppose that a list of L perhaps overlapping groups G 1 , ..., G L on the coefficients Z is given. That is, if we are learning a representation of X with K atoms, each G i ⊂ P ( 1 , ..., K , where P is the set of all subsets of its argument, is specified. We can generalize the LLoyd algorithm for K means or K flats to this setting. After initializing the dictionary W , we find the distance of each x in X to its projection P G i x onto the span of W G i for each i . Each x is associated to the i with the smallest distance

and we find the coefficients

Then we update W to be the minimum of the convex problem

Algorithm 1 SOMP [12] function Z = SOMP ( X,W,K ) Initialize: coefficients Z = 0 , residual R = X , active set Ω = ∅ . repeat j = argmax i ∑ s | W T i R s | Ω = Ω ⋃ j Z = ( W T Ω W Ω ) -1 W T Ω X R = X -WZ until K iterations end function

and repeat. Each of the subproblems either has an explicit solution or is convex, and so the energy decreases. When the training is finished, we define g to be the function that maps each point in x ∈ R d to the i minimizing the error of the projection of x onto the span of W G i .

We can also run the same sort of algorithm when in addition to each group G specifying a list of indices, it also specifies a cost for the use of the dictionary elements associated to each index. If we choose an l 2 cost for each of the coefficients, we still get explicit updates and the decrease of energy at each round.

Note that if the number of groups is very large, it may be too costly to find the best group for each x exaustively. However, we can make a greedy approximation by running a modified OMP. Here, supposing at iteration s of the OMPwehavean active set Ω , the available dictionary elements to add to Ω are the union of all groups containing Ω . It is not necessary to be able to enumerate all the groups to use this method, only to have a subroutine which given Ω ⊂ 1 , ..., k can return ⋃ Ω ⊂ G G . However, using this sort of greedy approximation removes the guarantee that the energy decreases at each iteration.

Learning the groups with simultaneous orthogonal matching pursuit

In the previous section the groups were specified in advance. If we want to learn the groups, we can add a step in the algorithm. Now instead of taking the list of groups as input, we instead input just the number K of dictionary elements and the number of coefficients allowed per x .

After associating to each x the group that best represents it, we can turn around and consider all the x associated to that group. Our task is then to choose a subset of the dictionary that best represents that group. A greedy approximation to this problem in the least squares sense is given by the Simultaneous Orthogonal Matching Pursuit algorithm (SOMP)[12]. This algorithm proceeds just as a standard OMP, but at each iteration, all the x associated to a given group have to choose the next dictionary element added to the group together. See algorithm 1.

Unfortunately, because neither OMP nor SOMP is guaranteed to find the optimal solution to the NP hard problems they address, the energy may not decrease at each iteration with this scheme; however, as usual, we have found that in practice these methods do usually lead to a decrease in the energy. As in K -means, it may happen that no group uses a dictionary element; in such a situation one can remove a dictionary element from one of the groups, find the residual, and replace the unused dictionary element by the principal component of the residual.

We note that the model presented here can be thought of as a greedy sparse coding version of a 'topic model'. The dictionary elements act as the words, the x as the documents, and the groups are the topics. The algorithm learns the topics and the dictionary simultaneously.

The main focus of this work will be choosing a g that can be computed rapidly and learning a dictionary that respects g . We will consider g to be a hash function on R d , and hash buckets will be the atomic units of the groups; that is, the groups will either be the hash buckets or will be glued together from the hash buckets. This can be considered a sort of geometric regularization of the sparse coding problem: the active set will be forced to remain constant on the region of R d corresponding to each hash bucket.

Once g is chosen, we will learn the dictionary (and perhaps groups) as above, but instead of allowing each x to choose the group that best represents it individually, the x in a hash bucket will need to choose the group that best represents them together on average. We will also try to approximate standard greedy dictionary learning; in this

Require: data X , number of dictionary elements K , number of active coefficients per data point q , number of iterations I , and if desired, g : R d ↦→ 1 , ..., M .

Experiments

• 1: Each x chooses a group in { 1 , ..., L } via (2) or by the modified OMP as in Section 2.1. If g is given, all the x in a hash bucket are forced to choose the same group.
• 2: Each group in { 1 , ..., L } chooses subset of { 1 , ..., K } using Z = SOMP ( X,W,k ) .
• 3: Update W , either via K -SVD, or a least squares solve.

until I iterations case, there will be one group for every hash bucket. As above, and as with K -means, it may happen that no spatial bucket uses a particular group; in that case we can just pick a bucket at random and use the output of SOMP on that bucket to regenerate the unused group.

Learning how to quantize R d is a much studied (but still not completely understood) problem. One common motivation is to build a data structure allowing nearest neighbors from a given data set to be quickly computed. Another common motivation is to use the buckets of the quantization as words to build bag of words feature representations. The relationship between vector quantization and sparse coding has studied before by many authors[]. In particular, K -means is simply l 0 sparse coding with only the coefficients 0 and 1 allowed, and only 1 nonzero per x 1 .

In this work we will use a 2 -means tree 2 with subdivisions along medians to define g . We start by taking the entire data set and running 2 -means, obtaining centers c 1 and c 2 . We take each data point x ∈ X and find the angle between x and c 1 -c 2 ; X is divided at the median. We then repeat on each of the pieces, continuing until each piece is within a given distance to its mean, or a set depth p , whichever comes first. We initialize the 2 -means with farthest insertion, as in [13]. Note that our

1 'shape gain coding' allows a non-binary coefficient

2 Although perhaps not exactly standard usage, we will call the data structure obtained from binary partitions of R d a hash

experience is that very few iterations are necessary, and really the farthest insertion is suffient; in fact cutting in random directions (with some additional tricks and randomizations) has been shown to lead to good partitions when the underlying data has a 'manifold' structure, see [14]. The number of buckets at the bottom of the tree is upper bounded by 2 p ; we will choose p small enough so that it is simple to store a lookup table with the indices into the dictionary for each bucket, as well as the decision vectors for each branch in the tree.

We also could use mappings of the form g ( x ) = s ( h ( Hx + b )) , where H is a p × d matrix, h is some sort of nonlinearity (e.g. tanh , or sin ), b is an offset, and s is a thresholding function []. These mappings require less storage and are somewhat simpler to compute for the same bit depth, but on the data sets we work on, they have the disadvantage that many of the buckets are often empty or have very few entries for reasonable p . While this can be remedied by simply gluing (nearly) empty buckets to nearby full buckets and updating the lookup table, we have found the trees to work better. Note also that unlike in nearest neighbor data structures, it is unnecessary for leaf nodes to keep track of spatially nearby leaves that are far away in the tree metric, because all we care about is which dictionary atoms are used at that node.

After building g and training the dictionary, in order compute the coefficients of a new data point x , we pass it though the tree, obtaining g ( x ) . We lookup g ( x ) in a table, and this gives an index of m columns Ω of W ; at this point we solve the linear system W Ω z -x to get the outputs. Alternatively, for each group, we can store ( W T Ω W Ω ) -1 (or some stable factorization), or ( W T Ω W Ω ) -1 W T Ω , and just do the requisite matrix multiplications

Discussion of related work

The idea of clustering the input space and then using a different dictionary for each cluster has appeared several times before. As mentioned in the introduction, a simple example is the K -flats algorithm, or other subspace clustering algorithms [15]. There, the subdictionaries serve the dual purpose of determining the clusters and also finding the coefficients for the data points associated to them. More recently this technique has been succesfully applied to object recognition by [6, 7]. In those works, the clus-

ters are determined by K -means (or a Gaussian mixture model); in the first, there is a different dictionary for each cluster, and the code is the size of the union of all the subdictionaries, but only the blocks corresponding to the centroids near the input are nonzero. In the second work, the dictionaries for each centroid are the same, but the code is still a concatenation of the codes associated to each centroid (and are set to zero if the input does not belong to that centroid). The current work differes from these in two ways. The first is the use of a fast method for clustering, and the second is the use of shared parts across the dictionaries, where the organization of the parts sharing has been learned from the data.

In [16] the authors construct a dictionary on the backbone of a hierarchical clustering with fast evaluation. They also use shared parts. However, in that work the part sharing is determined by the tree structure of the clustering, and not learned.

There is now a large literature on structured sparsity. Like this work, [11, 17] use a greedy approach for structured sparse coding based on OMP or CoSaMP. Unlike this work, they have provable recovery properties when the true coefficients respect the structure, and when the dictionaries satisify certain incoherence properites. On the other hand, those works do not attempt to learn the dictionary, and only discuss the forward problem of finding z from x and W .The works in [8, 9, 10] use an approach to structured sparsity that allows for convex optimization in z . In these works the coefficients are arranged into a predetermined set of groups, and the sparsity term penalizes the number of active groups, rather than the number of active elements; the dictionary is trained to fit the data. None of these works attempt to learn the group structure along with the dictionary

Finally we note that other works have explored the idea of accelerating sparse coding by training the dictionary along with an approximation method, e.g. [5, 18]. In the first, the approximation is via a single layer feed forward network, and in the second, via a multilayer feed forward network with a shrinkage nonlinearity. This work uses a tree and lookup table instead.

What do the groups look like?

To get a sense of what kind of groups learned from algorithm 2 look like, we train a dictionary on 500,000 8 × 8 image patches, and view the results. The image patches are drawn from the PASCAL dataset, and their means are removed. We train a dictionary with 256 elements and 512 groups; each group has 5 dictionary elements in it. We train using the batch method, with a K -SVD update for the dictionary.

After training, some of the dictionary elements are used by many groups, and others are used by only a few. The median number of groups using a given element is 6; 47 elements are in exactly 1 group, and 15 are in more than 30. In figure 1 we display the dictionary ordered by the number of groups containing each element; this number increases in each column and moving to the right. Unsurprisingly, 'popular' elements that belong to many groups are low frequency. In this figure we also show the groups containing a few chosen atoms.

Review of the image classification pipeline

Here we will review a standard pipeline for object recognition [1, 2], while giving details about our implementation, which streamlines certain components. It consists of the following parts: 1) Calculation of sift vectors at every location (sift grid) 2) Calculation of the feature vectors for every sift vector using the 'tree sparse coding' described above, 3) Spatial pyramidal max pooling 4) logistic regression or SVM classification. Care is taken to calculate each of these parts efficiently.

3.2.1 Sift grid

We run tests with two different implementations of dense sift. The first is matlab code by L. Lazebnik [1]. We also use a fast, approximate c++ version that we coded ourselves. The details are as follows:

The x and y derivatives. We convolve the image with two 5 × 5 filters that are the x and y derivatives of Gaussian. This results in the values of x and y derivatives

I y = dI/dy, I x = dI/dx of the image intensity at every location of the image.

Orientation histogram. This operation takes the two gradient values I y , I x at every location and smoothly bins them into histogram of eight orientations ( 0 , π/ 4 , . . . , 7 π/ 4 ) as follows. First we calculate the orientation angle φ = arctan( I y /I x ) + π (1 -sign ( x )) / 2 and magnitude m = √ I 2 x + I 2 y . Let φ h ( n ) = nπ/ 4 , n = 0 , . . . , 7 . The final set of values is v ( n ) = m ∗ cos( φ -φ h ( n )) 9 + where the x + = x if x > 0 and 0 otherwise. Most of these operations are computationally expensive and therefore we precompute these values. We bin the I y and I x values into 500 bins each and for every combination ( 500 2 values) we calculate v ( n ) , n = 0 , . . . , 7 . The bin range is chosen so that the values of I y and I x never fall outside the range of the binning so no checks are needed. After this computation we obtain 8 values at every location of the image.

Smooth subsampling We subsample the resulting features by two in each direction. Specifically let v n,y,x be the input value obtained from the previous step, where n is the feature number and y, x is the location. The output value will be u n,y,x = v n, 2 y, 2 x + v n, 2 y, 2 x +1 + v n, 2 y +1 , 2 x + v n, 2 y +1 , 2 x +1 . This is efficient since it only involves additions. Note that it results in output values that are essentially four times larger the input values at each location.

Smoothing We convolve each feature with [[1 , 1] , [1 , 1]] filter. This is calculated using u n,y,x = v n,y,x + v n,y,x +1 + v n,y +1 ,x + v n,y +1 ,x +1 again resulting in essentially four times larger output values then input values.

Combining and normalizing into sift vector Now we obtain 128 component sift vector from every location of the features maps from the previous step. At every location ( x, y ) (of the subsampled feature image) we first obtain 128 component vector by concatenating the 8 -component vectors at the following locations ( x +2 i, y + 2 j ) , i = 1 , 2 , 3 , 4 and j = 1 , 2 , 3 , 4 . Then we normalize this vector as follows. If the norm of the vector is smaller then the threshold t h = 1 we keep the vector. If it is larger we normalize it to have size t h . The result is placed into the appropriate location of the final m y × m x × 128 vector, where m y,x ≈ n y,x / 2 where n x,y are the dimen-

sions of the original image. The dimensions are slightly smaller due to boundary effects. This last operation (combining and normalizing) is the most expensive operation in the sift grid calculation and we took care to implement it efficiently. Note that in Lazebnik's (and Lowe's original) sift the smoothing is done over a larger neighborhood with inputs near the center weighted more then those further. This makes the output more smoothly varying under translations; in our case we used equal weighting over small neighborhoods for computational efficiency.

Hashed sparse coding.

We used the main procedure of this paper to calculate feature vector for each sift vector. Each such computation consisted essentially of depth= 16 multiplications of sift and tree decision vectors ( 16 × 128 computations) followed by multiplication of the sift vector by the appropriate pseudo-inverse matrix (typically 128 × 5 multiplications) resulting in total of approximately 128 × 21 multiplications. For 2048 dimensional feature vector this compares to 128 × (2048 + 4) multiplications that are needed for omp resulting in almost 100 -fold reduction. Our model was trained on 2 × 10 6 randomly selected sift vectors from Pascal 2011 dataset.

Spatial pyramidal pooling.

We used the same spatial pyramidal max pooling as in [YLan]. Since the feature vectors are in the sparse format the resulting computation is very efficient and negligible compared to either sift or tree sparse coding. The details are as follows. We need to calculate the maximum over the features in 1 × 1 , 2 × 2 and 4 × 4 regions of the feature vector obtained in the previous step. First we split this vector into 4 × 4 regions R x ′ ,y ′ . Let n f be number of features, typically 2048 , v f,x,y be the input feature vector and u f,x ′ ,y ′ , x ′ = 1 , 2 , 3 , 4 , y ′ = 1 , 2 , 3 , 4 be the 4 × 4 part of the final feature vector. We calculate u using the following.

This calculation is done by looping over all feature vectors and indices and filling the pooled feature vector so the number of computations is of the order of the total number of nonzero features. We can get 2 × 2 and 1 × 1 parts of the final feature vector analogously. However it is more efficient now to use the 4 × 4 vector obtained and pool it into 2 × 2 regions and then pool the result into 1 × 1 regions. The final output vector is concatenation of these vectors, resulting in n f × 21 vector.

Subsequently a logistic regression classifier is trained on the feature vectors using the liblinear package [19].

Each the following operations we implemented using a multicore processing: all steps of the sift, finding the group using tree, and multiplying by pseudo-inverses. In each of these steps separately the image/feature image was split in n cores parts and send to different core. The system was implemented in C++. Blas in the Accelerate framework was used in the tree sparse coding. We report the result on a macbook pro, with a 2.3 Ghz Intel Core i7 processor with 4 cores. The observed speedup compared to single core was about 3 .

We also test the run time of just the coding, compared with coding using OMP with the SPAMS package [20].

Accuracy on Caltech 101 and 15 scenes

We test the accuracy of the standard pipeline with the hashed dictionary and with standard l 0 sparse coding on two object recognition benchmarks, Caltech 101 [] and 15 scenes []. As mentioned before, for all data sets, we train the hashed dictionary on 2 × 10 6 randomly selected sift vectors from the Pascal 2011 dataset. Caltech 101 consists of 101 image categories and approximately 50 images per category; many classes have more training examples and we do the usual normalization of error by class size. We use 30 training examples per class. The 15 scenes database contains 15 categories and 4485 images, and between 200 to 400 images per category. We use 100 training images per class on this data set. For each data set, we run over 10 random splits and record the mean and standard deviation of the test error. We record the results in Tables 3.4 and 3.4. The first two columns of each table correspond to the hashed sparsed coding run with 5 or 10 nonzero entries on Lazebnik's sift. The next two columns

correspond to the 'real time' system, hashed sparse coding run on our approximate sift, and the last two columns correspond to OMP, trained and coded with SPAMS[20] on Lazebnik's sift. Each row corresponds to the number of atoms in the dictionary. As far as we know, state of the art with single features on grayscale images on Caltech 101 with 30 training examples per category is .773, in [7], and .898 for the 15 scenes, in [21]. Both of these methods use the same basic pipeline as this work, but with variations on the sparse coding; our method can be used in conjuction with their methods.

As has been observed by other authors, increasing the size of the dictionary only seems to increase the accuracy. Note that for our method, the only places that the size of the dictionary affects the computational cost is in training, where we use an SOMP, and in the final classification stage. The last component is small for these experiments, but if we wanted to use the system for detection at many locations at an image, it would start to be significant.

Running speed.

We tested the speed of the full pipeline from image to classification. We show results on images from the Berkeley dataset and Caltech 101. The Berkeley images are 321 × 481 , The Caltech 101 images were resized so that the largest size was at most 300, with the aspect ratio fixed. With 5 nonzero coefficients and depth 16 tree, we get the results in Table 3.4. The entire dataset of 9145 images in Caltech 101 was processed in 4 minutes and 48 seconds with 2048 features and in 5 minutes and 35 seconds with 8092 features. This corresponds to 31 . 75 fps and 27 . 3 fps respectively.

We also test the speed of just the sparse coding 3 . Coding 15000 sift vectors with a depth 16 tree and 5 nonzeros per x takes .034 seconds with one core, and .018 with four. In comparison, SPAMS with a dictionary of size 1024 costs .898 seconds using four cores. This is not exactly a fair test, as SPAMS must calculate a Cholesky decomposition of the Gram matrix of the dictionary when it runs, and this could be cached; however, simply multiplying the dictionary matrix by the data vectors takes .294

3 This test was done on a quad core intel i5 running 64 bit Linux, with 4 gigs of ram; both our code and SPAMS were run as a mex file through Matlab

seconds. As the size of the dictionary increases, this will increase, but our method will not get any slower.

Conclusion

In this paper we presented a fast approximate sparse coding algorithm and use it to build an accurate real time object recognition system. Our contributions can be summarized into four parts. 1) We describe a general method for learning the groups for greedy structured sparse coding using a generalization of LLoyd's algorithm and SOMP. 2) We use this method to design a fast approximation of greedy sparse coding that uses a tree structure for inference. 3) We give a fast approximate implementation of the SIFT descriptor. 4) These algorithms together allow as to build a real time object recognition system in the framework of [2]. It processes the entire Caltech 101 dataset in under 5 minutes (with images resized so that larger size is 300 pixels). As far as we know this is the first time that a fast implementation of this type of system has been put together with comparable accuracy.

We see many possible directions in the future both for improving the group sparse coding algorithm and applying our system to vision. We would like to learn the hash or tree, rather than build it before the dictionary training. We would like to train the system on larger datasets and work on real time object detection (as opposed to classification). At this speed the algorithm allows us to process around 2 million medium sized images ( 300 × 400 ) in 24 hours on a single computer. The object detection should also be feasible given that the expensive part calculation of features at different parts of the image from which detection is calculated - is fast.

Algorithm: SOMP [gisttr05]

$Z=SOMP(X,W,K)$
 
 Initialize: coefficients $Z=0$, residual $R=X$, active set $\Omega=\emptyset$.
 **repeat**
$j=\argmax_i \sum_s|W_i^TR_s|$
 $\Omega=\Omega \bigcup j$
 $Z=\left(W_\Omega^TW_\Omega\right)^-1W_\Omega^TX$
 $R = X-WZ$
 **until** $K$ iterations

Algorithm: Learning a dictionary and groups

Require: data $X$, number of dictionary elements $K$, number of active coefficients per data point $q$, number of iterations $I$, and if desired, $g:\R^d\mapsto \{1,...,M\}$.
 
 **repeat**
1: Each $x$ chooses a group in $\{1,...,L\}$ via or by the modified OMP as in Section . If $g$ is given, all the $x$ in a hash bucket are forced to choose the same group.
 2: Each group in $\{1,...,L\}$ chooses subset of $\{1,...,K\}$ using $Z=SOMP(X,W,k)$. 
 3: Update $W$, either via $K$-SVD, or a least squares solve. 
 **until** $I$ iterations

Table 3: Speeds of different parts of the system and different dictionary sizes on 321 × 481 pixel Berkeley dataset images and Caltech 101 images. The times are for single frame in seconds. Frame rates are the inverses and are in frames per second. The total time is the time to process the entire Caltech 101 datasets consisting of 9145 images (minutes:seconds). The Caltech 101 images were pre-resized so that largest side is 300 pixels. The last column is the recognition performance when trained on 30 training images per category. (The speeds vary probably due to disc access and are faster after one or more sweeps through the dataset).

Table: S3.T1: Caltech accuracies and standard deviations over 10 random splits. The first two columns of each table correspond to the hashed sparsed coding run with 5 or 10 nonzero entries, on Lazebnik’s sift. The next two columns correspond to the “real time” system, hashed sparse coding run on our approximate sift, and the last two columns correspond to OMP, trained and coded with SPAMS[20] on Lazebnik’s sift. Each row corresponds to the number of atoms in the dictionary.

	hashed m=5	hashed m=10	hashed m=5, R.T.	hashed m=10 R.T.	OMP m=5	OMP m=10
K=1024	.722 ± .011	.704 ± .010	.710 ± .007	.697 ± .010	.725 ± .008	.721 ± .010
K=2048	.735 ± .007	.731 ± .011	.723 ± .007	.716 ± .005	.747 ± .008	.738 ± .008
K=4096	.741 ± .011	.740 ± .006	.736 ± .005	.724 ± .004	.754 ± .008	.757 ± .010
K=8092	.751 ± .009		.739 ± .003

Table: S3.T2: 15 scenes accuracies and standard deviations over 10 random splits. The first two columns of each table correspond to the hashed sparsed coding run with 5 or 10 nonzero entries on Lazebnik’s sift. The next two columns correspond to the “real time” system, hashed sparse coding run on our approximate sift, and the last two columns correspond to OMP, trained and coded with SPAMS[20] on Lazebnik’s sift. Each row corresponds to the number of atoms in the dictionary.

	hashed m=5	hashed m=10	hashed m=5, R.T.	hashed m=10 R.T.	OMP m=5	OMP m=10
K=1024	.792 ± .006	.789 ± .004	.786 ± .004	.770 ± .007	.801 ± .006	.802 ± .004
K=2048	.807 ± .006	.800 ± .004	.796 ± .007	.788 ± .007	.814 ± .006	.813 ± .006
K=4096	.810 ± .007	.810 ± .004	.807 ± .003	.804 ± .004	.826 ± .007	.822 ± .007
K=8092	.811 ± .004		.815 ± .004

	321 × 481 pixel images	321 × 481 pixel images	321 × 481 pixel images	321 × 481 pixel images	Caltech 101 (on 4 cores)	Caltech 101 (on 4 cores)	Caltech 101 (on 4 cores)
	1 core (s)	4 cores (s)	1 core (fps)	4 cores (fps)	total time (m:s)	(fps)	performance
SIFT	0.039	0.017	25	59	-	-	-
SIFT+TreeSC+pyramid	0.143	0.045	7	22.5	-	-	-
full (1024)	0.145	0.0465	6.9	21	4:01	38	. 710 ± . 007
full (2048)	0.1473	0.050	6.8	20	4:45	32	. 723 ± . 007
full (4096)	0.1495	0.052	6.7	19	4:42	32	. 736 ± . 005
full (8092)	0.155	0.0565	6.4	18	5:35	27	. 739 ± . 003

References

[1] C. Schmid S. Lazebnik and J. Ponce, ``Beyond bags of features: Spatial pyramid matching for recognizing natural scene categories'', in CVPR'06, 2006.

[2] Jianchao Yang, Kai Yu, Yihong Gong, and Thomas Huang, ``Linear spatial pyramid matching using sparse coding for image classification'', in CVPR'09, 2009.

[3] B.A. Olshausen and D. Field, ``Emergence of simple-cell receptive field properties by learning a sparse code for natural images'', Nature, vol. 381, no. 6583, pp. 607--609, 1996.

[4] M. Aharon, M. Elad, and A. Bruckstein, ``K-SVD: An algorithm for designing overcomplete dictionaries for sparse representation'', IEEE Transactions on Signal Processing, vol. 54, no. 11, pp. 4311--4322, 2006.

[5] Koray Kavukcuoglu, Marc'Aurelio Ranzato, and Yann LeCun, ``Fast inference in sparse coding algorithms with applications to object recognition'', Tech. Rep. CBLL-TR-2008-12-01, Computational and Biological Learning Lab, Courant Institute, NYU, 2008.

[6] K. Yu J. Yang and T. Huang., ``Efficient highly overcomplete sparse coding using a mixture model.'', in European Conference on Computer Vision, 2010.

[7] Y. Boureau, N. La Roux, F. Bach, J. Ponce, and Y. LeCun, ``Ask the locals: multi-way local pooling for image recognition'', in International Conference on Computer Vision, 2011.

[8] R. Jenatton, J. Mairal, G. Obozinski, and F. Bach, ``Proximal methods for sparse hierarchical dictionary learning'', in International Conference on Machine Learning (ICML), 2010.

[9] Seyoung Kim and Eric P. Xing, ``Tree-guided group lasso for multi-task regression with structured sparsity'', in ICML, 2010, pp. 543--550.

[10] Laurent Jacob, Guillaume Obozinski, and Jean-Philippe Vert, ``Group lasso with overlap and graph lasso'', in Proceedings of the 26th Annual International Conference on Machine Learning, New York, NY, USA, 2009, ICML '09, pp. 433--440, ACM.

[11] Richard G. Baraniuk, Volkan Cevher, Marco F. Duarte, and Chinmay Hegde, ``Model-Based Compressive Sensing'', Dec 2009.

[12] Anna C. Gilbert, Martin J. Strauss, and Joel A. Tropp, ``Simultaneous Sparse Approximation via Greedy Pursuit'', IEEE Trans. Acoust. Speech Signal Process., vol. 5, pp. 721--724, 2005.

[13] R. Ostrovsky, Y. Rabani, L. Schulman, and C. Swamy, ``The effectiveness of lloyd-type methods for the k-means problem'', in FOCS 2006, 2006.

[14] S. Dasgupta and Y. Freund., ``Random projection trees and low dimensional manifolds.'', in STOC 2008, 2008.

[15] R. Vidal, ``Subspace clustering'', IEEE Signal Processing Magazine, vol. 28, pp. 52--68, 2011.

[16] W. Allard, G. Chen, and M. Maggioni, ``Multiscale geometric methods for data sets II: Geometric multi-resolution analysis'', to appear in Applied and Computational Harmonic Analysis.

[17] Junzhou Huang, Tong Zhang, and Dimitris N. Metaxas, ``Learning with structured sparsity'', in ICML, 2009, p. 53.

[18] K. Gregor and Y. LeCun, ``Learning fast approximations of sparse coding'', in International Conference on Machine Learning (ICML), 2010.

[19] Rong-En Fan, Kai-Wei Chang, Cho-Jui Hsieh, Xiang-Rui Wang, and Chih-Jen Lin, ``LIBLINEAR: A library for large linear classification'', Journal of Machine Learning Research, vol. 9, pp. 1871--1874, 2008.

[20] '', http://www.di.ens.fr/willow/SPAMS/.

[21] Liang-Tien Chia Shenghua Gao, Ivor Wai-Hung Tsang and Peilin Zhao, ``Local features are not lonely-laplacian sparse coding for image classification.'', in CVPR 2010, 2010.

---