Ala Aboudib

Vincent Gripon

Gilles Coppin

A MODEL OF BOTTOM-UP VISUAL ATTENTION USING CORTICAL MAGNIFICATION

Keywords: cortical magnification, bottom-up attention, saliency, multi-scale

The focus of visual attention has been argued to play a key role in object recognition. Many computational models of visual attention were proposed to estimate locations of eye fixations driven by bottom-up stimuli. Most of these models rely on pyramids consisting of multiple scaled versions of the visual scene. This design aims at capturing the fact that neural cells in higher visual areas tend to have larger receptive fields (RFs). On the other hand, very few models represent multi-scaling resulting from the eccentricity-dependent RF sizes within each visual layer, also known as the cortical magnification effect. In this paper, we demonstrate that using a cortical-magnification-like mechanism can lead to performant alternatives to pyramidal approaches in the context of attentional modeling. Moreover, we argue that introducing such a mechanism equips the proposed model with additional properties related to overt attention and distance-dependent saliency that are worth exploring.

INTRODUCTION

The study of visual attention and its computational modeling is an emerging field. Object recognition is one of its numerous potential domains of applications [START_REF] Borji | State-of-the-art in visual attention modeling[END_REF]. The seminal work of Treisman and Gelad in the Feature Integration Theory (FIT) [START_REF] Treisman | A featureintegration theory of attention[END_REF] investigated the key role of focusing attention in correctly associating visual properties related to the same object. Later on, Koch and Ullman [START_REF] Koch | Shifts in selective visual attention: towards the underlying neural circuitry[END_REF] introduced the first theoretical model suggesting a possible neuro-inspired bottom-up mechanism for guiding attention. The first working implementation of Koch and Ullman's model was later proposed by Itti and Koch in [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF] and became a landmark model for bottom-up attentional modelling.

Since then, many works aimed at dealing with the object recognition problem from an attentional perspective. Rybak [START_REF] Ia Rybak | A model of attentionguided visual perception and recognition[END_REF] proposed an algorithm that associates a scanpath drawn from several fixations with the contents of fixated regions to perform learning and recognition. The role of gist information and context in guiding visual attention was elaborately studied by Torralba [START_REF] Torralba | Contextual guidance of eye movements and attention in real-world scenes: the role of global features in object search[END_REF], leading to a better understanding of the different kinds of top-down influences on attention. Later on, Walther and Koch proposed a unified framework for combining attention and object recognition [START_REF] Dirk | Attention in hierarchical models of object recognition[END_REF].

The main contribution of this paper is to demonstrate that cortical-magnification-like mechanisms can be applied as an alternative to pyramidal multi-scale processing traditionally used in computational models of visual attention. We introduce a new model and we compare it to the renowned Itti and Koch model originally proposed in [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF]. We obtain similar or even better performance while using a much smaller number of feature maps. After that, we explore some interesting properties the magnification mechanism adds to the proposed model.

The rest of this paper is organized as follows: Section 2 introduces some related work. In section 3 we explain the technical details of the proposed model. We introduce and discuss simulation results in section 4. Section 5 is a conclusion.

RELATED WORK

Image pyramids are widely used tools in models of bottom-up visual attention. For instance, Itti and Koch in their famous work [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF] used dyadic Gaussian pyramids to generate several spatial scales of each input image. Feature maps are then created by applying a cross-scale subtraction operator on these images.

A different method was proposed by Rybak in [START_REF] Ia Rybak | A model of attentionguided visual perception and recognition[END_REF]. He modelled multiple spatial scales by creating a single retinal image centered at each fixation point. Each such image is composed of three concentric regions with increased blurring as they go far from the fixation center. Advani in his his MR-AIM model [START_REF] Advani | A multi-resolution saliency framework to drive foveation[END_REF] creates multi-resolution images by concatenating patches extracted from a Gaussian pyramid. Multiresolution processing was also used in many other models such as [START_REF] Marfil | Combining segmentation and attention: a new foveal attention model[END_REF] and by using wavelet transforms in [START_REF] Li | A saliency model based on wavelet transform and visual attention[END_REF].

More recently, the computational modeling of cortical magnification began to emerge in a more formal way. One example is the work of Freeman and Simoncelli on Metamers in the ventral stream [START_REF] Freeman | Metamers of the ventral stream[END_REF]. Another example is the work of Poggio [START_REF] Poggio | Computational role of eccentricity dependent cortical magnification[END_REF] on the computational role of cortical magnification as a sampling extension of the "Magic theory" [START_REF] Anselmi | Magic materials: a theory of deep hierarchical architectures for learning sensory representations[END_REF] and its role in scale and shift invariance for pattern recognition [START_REF] Isik | A hierarchical model of peripheral vision[END_REF].

METHODOLOGY

The proposed model is based on the original work of Itti and Koch in [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF]. Figure 1 shows the basic architecture of the proposed model. Input is provided to the model as an RGB image M . It is then embedded into a square black background to form the RGB scene E. We do such embedding in order to force the input to be square-shaped for an easier manipulation. Then, a grey-scale version of the scene E r is created as a simple average of the R,G and B channels of E.

A visual angle Θ M is associated with the image M . This simulates the fact that viewing a given image from different distances changes the visual angles occupied by that image and thus might change the attentional behavior of the visual system. Another visual angle Θ v is assigned to a central region of E that represents the fovea. Thus, the "further" the image is from the model's "eye", the smaller is the visual angle Θ M occupied by the image M and the larger is the area that falls within the fovea. This is meant to simulate the case when an image is seen by a human subject from different distances.

We do not use Gaussian pyramids to generate multiple scales of the input scene as in [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF]. We replace this step by generating feature maps using kernels with eccentricitydependent sizes representing RFs. Eccentricity dependency is observed across several layers of the ventral stream including LGN, V1, V2, V4 and beyond [15][16]. This allows us to reduce the number of required feature maps from 42 in the model of Itti and Koch [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF] to 9 maps in the proposed model while gaining in performance on predicting fixations (cf. Section 4).

Designing the kernel bank F

In order to create feature maps, a bank of spatial filters F is applied to the scene E and to its gray-scale version E r . Filter kernels in F are designed to emulate the eccentricitydependency of RFs in LGN and V1. This dependency is behind what is known as the cortical magnification effect. Therefore, kernels in the center of our model's "visual field" MVF should be more numerous and smaller in size and become sparser with a linear increase in their diameters as they approach the periphery [START_REF] Gattass | Visual topography of v2 in the macaque[END_REF] [START_REF] Gattass | Visuotopic organization and extent of v3 and v4 of the macaque[END_REF]. To approximate this effect, we determine the locations of kernel centers and their sizes using a geometry very similar to the foveal Cartesian geometry proposed by [START_REF] Martínez | A new foveal cartesian geometry approach used for object tracking[END_REF] as explained below. A square 2D map with the same height and width as the scene E is created as in the illustrative Figure 2(left). A central square region is chosen to represent the fovea containing a fixed number of kernel centers. Then, a number of concentric square rings each made of equally-spaced kernel centers on each side is created for the periphery. The side of each ring contains two more centers than its direct preceding inner ring. This allows these concentric rings to fit into a smaller square map which has as many pixels as there are kernels in F .(cf. Figure 2(right)). We refer to the latter as a "magnified map" and a subscript g is associated to its name. When a given kernel is applied to the scene by an inner product, the resulting value is assigned to the corresponding pixel in the magnified map.

All kernels whose centers are within the the model's fovea have the same radius. Similarly, all periphery kernels whose centers belong to the same square have the same radius. These radii are calculated from the linear equation r = ax where x is the eccentricity (in visual angles) of the corner pixels of the square to which a kernel belongs and a is the inverse magnification factor M -1 .

The eccentricity of corner pixels of the fovea are determined from the visual angle Θ v assigned to the fovea. Eccentricities of corner pixels of periphery rings are determined in such a way that an overlap p occurs between kernels of F along the diagonal.

Magnified feature maps

Nine feature maps are directly calculated; one map for the intensity modality, 4 maps for the color modality and 4 maps for the orientation modality.

Only one magnified intensity feature map I g is created from the gray-scale scene E r . A filter bank F DoG is used with L2-normalized Difference of Gaussian (DoG) kernels with zero mean to approximate the center-surround functionality of RFs in LGN. Then we apply a special convolution operator r between the gray-scale scene E r and the filter bank F DoG :

I g = E r r F DoG (1)
The operator r consists of a rectified normalized inner product between each kernel in F DoG and the patch it spans in E r . The scalar value resulting from the inner product then occupies one pixel in the magnified map I g as explained in Section 3.1.

Four magnified orientation feature maps O g (θ) with θ ∈ {0 • , 45 • , 90 • , 135 • } are created. Then the same convolution operator is applied:

O g (θ) = E r r F O(θ) (2)
where F O(θ) is a filter bank F with L2-normalized Gabor kernels [START_REF] Gabor | Theory of communication. part 1: The analysis of information[END_REF] with orientation θ and zero mean. Each Gabor filter f O(θ) ∈ F O(θ) is characterized by its orientation θ, effective width σ, wavelength λ and its aspect ratio γ. Finally, we use F DoG and the RGB scene E to create four magnified feature maps for the color-opponent channels RG g , GR g , BY g , YB g as follows:

RG g = E rg F DoG (3) GR g = E gr F DoG (4) BY g = E by F DoG (5) YB g = E yb F DoG (6)
Each of the convolution operators rg , gr , by and yb consists of a half-wave rectified (negative values are set to zero) normalized inner product between each kernel in F DoG and a corresponding 2D patch t. In the case of rg , the patch t is a concatenation between the area of the R component of E spanned by the center of the kernel and the area of the G component corresponding to the surround. The inverse holds for gr . Similarly, when applying the operator by , the patch t becomes a concatenation between the center area of the B component with the surround area of the yellow component Y (the average of the red and green components). The inverse hold for yb . This transformation is inspired by the DKL color space.

It should be noted here that kernel filters that overlap with the image M borders embedded in the scene E or E r could sometimes give high saliencies. This is most likely a result of the saliency of the image relative to the black background in which it is embedded. So all such values are set to zeros in all feature maps.

Magnified conspicuity maps and the saliency map

Three magnified conspicuity maps Īg , Ōg and Cg are created for intensity, orientation and color, respectively. The process is as follows:

Īg = N (I g ) (7
)

Ōg = N   θ∈{0 • ,45 • ,90 • ,135 • } N (O g (θ))   (8) Cg = N [N (|RG g -GR g |) + N (|BY g -YB g |)] (9)
where N (.) is the normalization operator designed by Itti and Koch in [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF]. The global magnified saliency map is then calculated as a linear combination of conspicuity maps:

S g = 1 3 ( Īg + Ōg + Cg ) (10)
Then a simple inverse mapping of pixels in S g into a map with the same size as the scene E is done. After this, the black background is removed to get a standard saliency map S that has the same size as the input image M . Notice that S is a very sparse map composed of a central fovea surrounded by concentric square rings. In order to be able to compare to other models, a continuous saliency map is created by convolving a Gaussian kernel on the most salient locations extracted from S using an alternating Winner-Take-All (WTA) and Inhibition of Return mechanism. See Figure 3.

RESULTS AND DISCUSSION

We used the Judd Benchmark [START_REF] Judd | A benchmark of computational models of saliency to predict human fixations[END_REF] for partially optimizing some parameters such as the number of fixations extracted, the width of the Gaussian kernel convolved on fixation points to create the final saliency map and the value of Θ M . Then, we ran the model on the MIT Saliency Benchmark [START_REF] Bylinskii | Mit saliency benchmark[END_REF] which contains 300 natural color images. We set the visual angle Θ M = 50 • for each input image M and Θ v = 1 • for the fovea which is close to its actual size in the human eye. The number of kernels in a bank F corresponding to the fovea is set to 41x41 units which is close to the number of units in the foveola of V1 estimated in [START_REF] Poggio | Computational role of eccentricity dependent cortical magnification[END_REF]. The overlap parameter p is set to 0.8. The slope of the inverse magnification curve presented in Section 3.1 is set to a = 0.16 as suggested by Gattass in [START_REF] Gattass | Visuotopic organization and extent of v3 and v4 of the macaque[END_REF] for V1. Parameters of Gabor kernels are inspired by [START_REF] Serre | Robust object recognition with cortex-like mechanisms[END_REF]; we set γ = 0.3, σ/λ = 0.8 and d/σ = 2.5 where d is the diameter of a given kernel. Parameters of center-surround DoG kernels were adapted from [START_REF] Robert W Rodieck | Quantitative analysis of cat retinal ganglion cell response to visual stimuli[END_REF]; we set d/σ c = 11.8, σ s /σ c = 3 and g s /g c = 0.8 where σ c and σ s are the standard deviations of the center and surround Gaussians of a given DoG kernel, respectively. While g c and g s are the strengths of the center and surround Gaussians, respectively. Table 1 shows the scores of our models according to several metrics used by the benchmark and how they are compared to the Itti and Koch model [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF]. Our results and comparisons with many more models are also available on the MIT Saliency Benchmark website http://saliency.mit.edu.

For all metrics, the proposed model performes much better than IttiKoch1 model which is the original implementation of [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF] by its authors. On the other hand, the IttiKoch2 implementation of the same model by Jonathan Harel as a part of his GBVS toolbox gives a roughly similar performance; It has a similar performance for the similarity metric and higher AUC (Area Under the ROC curve) metrics. However, the proposed model has a better Correlation-Coefficient (CC), Normalized Scanpath Saliency (NSS) and Earth Mover's Distance (EMD).

Metric

Our It is worth noting that the IttiKoch2 implementation and many other models are optimized for blur and center-bias. The proposed model, however, has no explicit center-bias applied and only a minor blur optimization. We hypothesize that an effect similar to center-bias arise naturally in our model due the magnification factor. This is more biologically plausible than explicitly applying a Gaussian mask to account for center-bias.

CONCLUSION

In this paper, we introduced a new model of bottom-up visual attention based on Itti and Koch work in [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF]. We demonstrated that using cortical magnification could be an alternative to pyramidal approaches to multi-scale image processing for modeling attention. We also showed that this allows us to boost performance while using a much smaller number of feature maps which can automatically account for center-bias phenomena.

The proposed model has also several other properties discussed below. They are worth investigating in future works.

Associating a visual angle Θ M to input images allows us to study the effect of the distance from which an image is displayed on saliency. This is an interesting effect that is not often considered by saliency models or benchmarks. We suggest that distance-dependent saliency should be evaluated by future benchmarks.

The proposed model is inherently adapted to performing overt attention where the position of the fovea moves across the scene. This is not directly the case for most attentional models. Given the magnification property, different saliency maps of the same scene are generated each time the fovea moves. This could have important implications on rethinking the role of the inhibition of return.

Finally, kernels with different sizes provide different information. Larger kernels provide more global information about the scene. By exploiting this property and using kernels as inputs to associative memories, one can imagine a unified framework where attention, gist information and object recognition techniques can all interact together.

  This work was supported by the European Research Council under the European Union's Seventh Framework Program (FP7/2007-2013) / ERC grant agreement n[Pleaseinsert"PrerenderUnicode-˝intopreamble] 290901.

Fig. 1 :

 1 Fig. 1: Model architecture.

Fig. 2 :

 2 Fig. 2: Placement of kernel centers (left) and mapping to a magnified map (right). Figure adapted from [18].

Fig. 3 :

 3 Fig. 3: Fixations generated by the proposed model and the corresponding saliency maps: (a) Image M (from the MIT Saliency Benchmark) with the first 50 fixations (b) The magnified saliency map S g (c) A cropped region of the saliency map S after inverse mapping from S g showing the fovea and some periphery square rings (d) The continuous saliency map made with a Gaussian convolved on fixation locations.

Table 1 :

 1 Comparison between the performance of the proposed model and Itti and Kock model on the MIT Saliency Benchmark.

		model IttiKoch1 IttiKoch2
	Similarity	0.44	0.20	0.44
	AUC (Judd)	0.74	0.60	0.75
	AUC (Borji)	0.72	0.54	0.74
	Shuffled AUC	0.58	0.53	0.63
	CC	0.39	0.14	0.37
	NSS	0.99	0.43	0.97
	EMD	4.24	5.17	4.26


Member, IEEE Huicheng Zheng

Grégoire Lefebvre

Christophe Laurent

Fast-Learning Adaptive-Subspace Self-Organizing Map: An Application to Saliency-Based Invariant Image Feature Construction

Keywords: Adaptive-Subspace Self-Organizing Map, feature construction, adult image filtering, image classification

The Adaptive-Subspace Self-Organizing Map (AS-SOM) is useful for invariant feature generation and visualization. However, the learning procedure of the ASSOM is slow. In this paper, two fast implementations of the ASSOM are proposed to boost ASSOM learning based on insightful discussions of the basis rotation operator of ASSOM. We reveal the objective function approximately maximized by the classical rotation operator. We then explore a sequence of two schemes to apply the proposed ASSOM implementations to saliency-based invariant feature construction for image classification. In the first scheme, a cumulative activity map computed from a single ASSOM is used as descriptor of the input image. In the second scheme, we use one ASSOM for each image category and a joint cumulative activity map is calculated as the descriptor. Both schemes are evaluated on a subset of the Corel photo base with 10 classes. The multi-ASSOM scheme is favored. It is also applied to adult image filtering and shows promising results.

I. INTRODUCTION

T HE Adaptive-Subspace Self-Organizing Map (ASSOM) proposed by Kohonen [START_REF] Kohonen | The Adaptive-Subspace SOM (ASSOM) and its use for the implementation of invariant feature detection[END_REF], [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF] is basically a combination of the SOM [START_REF] Kohonen | Self-Organizing Maps[END_REF] and the subspace method. By setting filters to correspond to basis vectors that span pattern subspaces, some transformation groups can be taken into account automatically. The ASSOM is an alternative to the standard Principal Component Analysis (PCA) method of feature generation. An earlier neural approach for PCA can be found in [START_REF] Oja | Principal components, minor components, and linear neural networks[END_REF]. Equivalence between probabilistic PCA and a typical subspace method for Gaussian density estimation has been established [START_REF] Wang | Links between PPCA and subspace methods for complete Gaussian density estimation[END_REF]. The ASSOM can generate spatially ordered feature filters thanks to interactions among processing modules [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF]. The input to an ASSOM array is typically an episode, i.e. a sequence of pattern vectors supposed to approximately span certain subspace. Typical examples of episodes include sequences of temporally consecutive speech signals or of image patches subject to transformations. By learning the episode as a whole, the Manuscript received June 20, 2006. This work was carried out during the tenure of a MUSCLE Internal fellowship (http://www.muscle-noe.org).

H. Zheng is with Media Processing and Communication Lab, Department of Electronics and Communication Engineering, Sun Yat-sen University, 510275 Guangzhou, China (email: zhenghch@mail.sysu.edu.cn) G. Lefebvre and C. Laurent are with France Telecom R&D -TECH/ IRIS/CIM, 4 Rue du Clos Courtel, 35512 Cesson Sévigné Cedex, France (email: {gregoire.lefebvre, christophe2.laurent}@orange-ftgroup.com)

Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org.

ASSOM is able to capture the transformation coded therein. The simulation results in [START_REF] Kohonen | The Adaptive-Subspace SOM (ASSOM) and its use for the implementation of invariant feature detection[END_REF] and [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF] have demonstrated that the ASSOM can induce ordered filter banks to account for translation, rotation and scaling. The relationship between the neurons in the ASSOM and their biological counterparts are reported in [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF]. The ASSOM has been successfully applied to speech processing [START_REF] Hase | Speech signal processing using Adaptive Subspace SOM (ASSOM)[END_REF], texture segmentation [START_REF] Ruiz-Del-Solar | TEXSOM: Texture segmentation using Self-Organizing Maps[END_REF], image retrieval [START_REF] De Ridder | The Adaptive Subspace Map for image description and image database retrieval[END_REF] and image classification [START_REF] De Ridder | The Adaptive Subspace Map for image description and image database retrieval[END_REF], [START_REF] Zhang | Handwritten digit recognition by Adaptive-Subspace Self-Organizing Map (ASSOM)[END_REF]. A supervised ASSOM was proposed by Ruiz-del-Solar in [START_REF] Ruiz-Del-Solar | TEXSOM: Texture segmentation using Self-Organizing Maps[END_REF].

The traditional learning procedure of the ASSOM involves computations related to a rotation operator matrix, which not only is memory demanding, but also has a computational load quadratic to the input dimension, i.e. the dimension of input vectors. Therefore, this algorithm in its original form is costly for practical applications, especially for image processing, where input patterns are often high dimensional. In order to reduce the learning time, the Adaptive Subspace Map (ASM) proposed by De Ridder et al. [START_REF] De Ridder | The Adaptive Subspace Map for image description and image database retrieval[END_REF] drops topological ordering and performs a batch-mode updating of the subspaces with PCA. However, without topological ordering, it is no longer ASSOM. López-Rubio et al. [START_REF] López-Rubio | A Principal Components Analysis Self-Organizing Map[END_REF] proposed the PCASOM by combining PCA with ASSOM, which runs about twice faster than the basic ASSOM under similar classification performance. López-Rubio et al. [START_REF] López-Rubio | New learning rules for the ASSOM network[END_REF] proposed two new learning rules of the ASSOM based on a gradient-based approach or on the Levenberg-Marquardt method. The new rules converge faster than the traditional ASSOM. However, the gradient-based approach showed obvious oscillations in the learning curve. The Levenberg-Marquardt method, on the other hand, has a computational load cubic with the input dimension, which excludes its use for high-dimensional inputs. McGlinchey et al. [START_REF] Mcglinchey | Fast formation of invariant feature maps[END_REF] replaced the traditional basis vector updating formula with one proposed by Oja [START_REF] Oja | Neural networks, principal components, and subspaces[END_REF], where the computational load is only linear to the input dimension, but quadratic to the subspace dimension. The above-mentioned methods can beat the ASSOM realized in the traditional way in terms of learning speed, but not the fast implementations of the ASSOM to be discovered in this paper.

In this paper, we shall show that in the ASSOM learning, the increment of each basis vector is a scaling of the component vectors of the input episode, which leads to a computational load linear to both the input dimension and the subspace dimension. We reveal that the rotation operator proposed in the original ASSOM approximately maximizes another objective function. The related in-depth analysis leads to a batch-mode 0000-0000/00$00.00 c 2008 IEEE updating of the basis vectors, where the increment of each basis vector is a linear combination of component vectors in the episode. This modified operator further accelerates the learning procedure.

In this paper, we are also interested in applying the fast implementations of ASSOM to saliency-based invariant feature construction for image classification. Saliency-based approaches are driven by psychovisual works [START_REF] Gordon | Theories of Visual Perception[END_REF], [START_REF] Marr | [END_REF], which have shown that the sensitivity of the human visual system (HVS) is not uniformly distributed across the image content. Such local methods can focus on different "concepts" in the image. Research has become very active on saliency-based approaches [START_REF] Schmid | Local grayvalue invariants for image retrieval[END_REF], [START_REF] Lowe | Object recognition from local scale-invariant features[END_REF], [START_REF] Schmid | Evaluation of interest point detectors[END_REF], [START_REF] Smeulders | Content-based image retrieval at the end of the early years[END_REF], [START_REF] Sebe | Salient points for content-based retrieval[END_REF], [START_REF] Mikolajczyk | A performance evaluation of local descriptors[END_REF]. Salient features can be represented by salient regions resulting from an image segmentation step [START_REF] Carson | Blobworld: Image segmentation using expectation-maximization and its application to image querying[END_REF], [START_REF] Wang | SIMPLIcity: Semantics-sensitive integrated matching for picture libraries[END_REF], by non-connected image zones resulting from the construction of saliency maps [START_REF] Itti | A model of saliency-based visual attention for rapid scene analysis[END_REF], by edges [START_REF] Gevers | Pictoseek: combining color and shape invariant features for image retrieval[END_REF], or by special points [START_REF] Harris | A combined corner and edge detector[END_REF], [START_REF] Mikolajczyk | An affine invariant interest point detector[END_REF]. In this paper, we are mainly interested in salient points, which provide the most compact representation of the image content by limiting the correlation and redundancy [START_REF] Laurent | In-depth analysis and evaluation of saliency-based color image indexing methods using wavelet salient features[END_REF]. The corner detector proposed by Harris and Stephens [START_REF] Harris | A combined corner and edge detector[END_REF] is one of the most used and the most known point detectors. It is based on the computation of a local auto-correlation function at each pixel location. Eigenvalues of the Hessian matrix of this function are used as indicator of presence of corners. In [START_REF] Laurent | In-depth analysis and evaluation of saliency-based color image indexing methods using wavelet salient features[END_REF], Laurent et al. proposed a wavelet-based salient point detector, which aims to address the following two observations: 1) contours are more perceptually important than other point related features such as corners; 2) salient points detected by a corner detector may be gathered in small image regions in the case of textured or noisy images, resulting in a very local image description. This wavelet-based detector reaches photometric invariance by incorporating a color invariance method proposed in [START_REF] Geusebroek | Color invariance[END_REF].

Various descriptors have been studied for saliency-based approaches [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF], [START_REF] Mikolajczyk | A performance evaluation of local descriptors[END_REF]. Lowe proposed the Scale Invariant Feature Transform (SIFT) for transform of image data into scaleinvariant local features [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF]. The SIFT descriptor was used in [START_REF] Csurka | Visual categorization with bags of keypoints[END_REF] for image classification, where images are represented by "bags of keypoints" through clustering the SIFT descriptors and the Support Vector Machine (SVM) is finally implemented for classification. Ideally, salient features should be robust to geometric transformations, slight changes of viewpoint and variations of imaging conditions, which is hard for handcrafted feature extractors. The ASSOM is a suitable tool to deal with transformations by learning pattern subspaces. The previous applications of the ASSOM to invariant feature extraction [START_REF] Ruiz-Del-Solar | TEXSOM: Texture segmentation using Self-Organizing Maps[END_REF], [START_REF] De Ridder | The Adaptive Subspace Map for image description and image database retrieval[END_REF], [START_REF] Zhang | Handwritten digit recognition by Adaptive-Subspace Self-Organizing Map (ASSOM)[END_REF] were not based on saliency approaches. In this paper, we explore a sequence of two schemes to apply the ASSOM to invariant feature construction for image classification under a saliency framework. In the first scheme, a single ASSOM is trained on all image patches extracted at the salient points from the training set by using the wavelet-based detector proposed in [START_REF] Laurent | In-depth analysis and evaluation of saliency-based color image indexing methods using wavelet salient features[END_REF]. An image is represented by a cumulative activity map through projecting the patches extracted from the image on this ASSOM. In the second scheme, one ASSOM per class is trained and an image is represented by a joint cumulative activity map calculated from these ASSOMs. For both schemes, the SVM will then be implemented for classification. It may be interesting that in [START_REF] Csurka | Visual categorization with bags of keypoints[END_REF], feature invariance is obtained by using SIFT descriptors before feature clustering, whereas in this paper, feature invariance is achieved at the clustering stage, by using the ASSOM.

c i (a) 1 b 2 b H b Q x 1 T b x 2 ˆL x 2 T b x H b x T (b)
This paper is organized as follows. In Section II, we discuss the traditional ASSOM learning algorithm and present the alternative fast-learning implementations. The performance in terms of learning speed will be demonstrated by experiments. In Section III, we will apply the ASSOM to saliency-based invariant feature construction for image classification. The experimental results will be presented in Section IV. This paper will be concluded by Section V, which summarizes main points in the paper and give some perspectives.

II. BASIS LEARNING RULES OF THE ASSOM

A. Basic ASSOM Learning

An ASSOM is composed of an array of modules with each one being realized by a two-layer neural network [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF], as shown in Fig. 1. In the two dimensional case, a lattice with a rectangular or hexagonal neighborhood is often used as the layout of modules, as in the SOM. Fig. 1(a) shows an ASSOM lattice with a hexagonal neighborhood. The architecture of a neural module of ASSOM is illustrated in Fig. 1(b). Let L be a subspace spanned by H orthonormal basis vectors {b 1 , b 2 , . . . , b H } and x be an input vector. The neurons in the first layer compute the orthogonal projections x T b h of x on the individual basis vectors b h , for h = 1, . . . , H. The only quadratic neuron of the second layer sums up the squared outputs of the first-layer neurons. The output of the module is then xL 2 , with xL being the orthogonal projection of x on L. It can be regarded as a measure of the matching between the input vector x and the subspace L. For an input episode X = {x(s), s ∈ S}, where S is the index set of vectors in the episode, Kohonen proposed to use the energy s∈S xL (s) 2 as the measure of matching between X and L [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF].

Modules in the ASSOM compete with each other based on their output energies. Once the winner of the modules is determined, the winning module and its neighbors update their subspaces to better represent the input subspace. A neighborhood function h (i) c is defined on the lattice of modules, where c indicates the index of the winning module and i the index of an arbitrary module. h

(i)
c is a function of learning step and the area of neighborhood shrinks with the learning step.

The classical Kohonen's ASSOM learning algorithm can be summarized as follows. At the learning step t, 1) Input the episode x(s), s ∈ S. Locate the winning module indexed by c = arg max i∈I s∈S xLi (s) 2 , where I is the set of indices of modules in the ASSOM. 2) For each module i in the neighborhood of c, including c itself, and for each input vector x(s), update the basis vectors b

(i)
h , according to the following procedure: a) Rotate each basis vector according to:

b (i) h = P (i) c (s, t)b ′ (i) h , (1)
where b

(i)
h is the new basis vector after rotation and b ′ (i) h the old one. P (i) c (s, t) is a rotation operator matrix defined by:

P (i) c (s, t) = I+λ(t)h (i) c (t) x(s)x T (s) xLi (s) x(s) . (2)
λ(t) is a learning-rate factor that diminishes with t.

For the algorithm to converge, λ(t) should satisfy

t λ(t) = ∞ and t λ 2 (t) < ∞ [32]. b) Dissipate the components b (i)
hj of the basis vectors b (i) h to improve the stability of the results [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF]: b

∼(i) hj = sgn(b (i) hj ) max(0, |b (i)
hj | -ε), where ε is a small positive value. c) Orthonormalize the basis vectors in module i. A naive implementation of (1) and ( 2) requires a matrix multiplication which needs not only a large amount of memory, but also a computational load quadratic to the input dimension. It would be costly for practical applications.

B. Insight on the Basis Vector Rotation

In this section, we propose an alternative implementation of the basis updating rule in the ASSOM learning. In the first place we propose to reformulate (1) and [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF]. The term b ′ (i) h in (1) can be distributed to the right side of (2), leading to:

b (i) h = b ′ (i) h + ∆b (i) h , (3)
where

∆b (i) h = λ(t)h (i) c (t) x(s)x T (s)b ′ (i) h xLi (s) x(s) . (4)
x T (s)b

′ (i) h
is in fact a scalar value. The equation becomes:

∆b (i) h = α (i) c,h (s, t)x(s) . (5)
Here α

(i)
c,h (s, t) is a scalar value defined by:

α (i) c,h (s, t) = λ(t)h (i) c (t) x T (s)b ′ (i) h xLi (s) x(s) . ( 6
)
This shows that the increment ∆b

(i)
h is in fact a scaling of the component vector x(s), as illustrated in Fig. 2, which seems to have been ignored by many practitioners. Careful examination of (5) would reveal similarity of this formula with a recursive PCA suggested in [START_REF] Oja | Subspace Methods of Pattern Recognition[END_REF]. The main difference here is that the gain of stochastic approximation is modulated by a neighborhood function dependent on module competition. Note that in (6), x T (s)b

) (i h b ) ( ' i h b ) ( ) , ( ) ( , ) ( s t s i h c i h x b α = ∆ ) (s x
′ (i)
h is the projection of the component vector on the basis vectors represented by the neurons of the first layer, which we have already when computing the projection xLi (s) (cf. Fig. 1). If we calculate the scaling factor α (i) c,h (s, t) first, and then scale the component vector x(s) with this factor, the computations associated with the basis vector updating will be dramatically reduced. This implementation will be referred to as FL-ASSOM for fast-learning ASSOM. It is completely equivalent to the basic ASSOM in terms of generating topologically ordered invariant-feature filters. Now we compare the computational loads of the basis vector updating in the basic ASSOM and in the FL-ASSOM. Let N be the input dimension. It is not hard to verify that a naive implementation of the updating rule defined by the matrix multiplications in ( 1) and ( 2) would need about HN 2 + N 2 scalar multiplications and about the same number of scalar additions. So the computational load is approximately O(HN 2 ), i.e. quadratic to the input dimension and linear to the subspace dimension. The replacement proposed by McGlinchey et al. [START_REF] Mcglinchey | Fast formation of invariant feature maps[END_REF] leads to a computational load of O(H 2 N ), i.e. linear to the input dimension but quadratic to the subspace dimension. Now with the proposed updating rule, the computations of xLi (s) and x(s) in ( 6) need about HN + 2N multiplications, and α (i) c,h (s, t)x(s) in ( 5) about HN multiplications. In all ( 5)and ( 6) need about 2HN + 2N multiplications. Similarly, the number of additions can be shown to be about 2HN + 2N . So with ( 5)and ( 6), the computational load is approximately O(HN ), i.e. linear to both the input dimension and the subspace dimension. The methods proposed in [START_REF] De Ridder | The Adaptive Subspace Map for image description and image database retrieval[END_REF], [START_REF] López-Rubio | A Principal Components Analysis Self-Organizing Map[END_REF], [START_REF] López-Rubio | New learning rules for the ASSOM network[END_REF], [START_REF] Mcglinchey | Fast formation of invariant feature maps[END_REF] can beat the classical implementation of the ASSOM, but not the implementation proposed here.

C. Discussion on the Operator

P (i) c (s, t)
The objective function of the ASSOM is defined as:

E = i h (i) c s∈S xLi (s) 2 x(s) 2 P (X)dX (7)
where P (X) is the distribution function of the random episode X. xLi (s) = x(s)-xLi (s) is the residual of x(s) after projection on L i . In the following we review briefly the main steps in the derivation of

P (i)
c (s, t) from ( 7) and prove that P (i) c (s, t) approximately maximizes another objective function.

By using the Robbins-Monro stochastic approximation [START_REF] Robbins | A stochastic approximation method[END_REF], a sample objective function has been aimed at:

E s (t) = i h (i) c (t) s∈S xLi (s) 2 x(s) 2 (8)
Kohonen showed in [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF] that

∂E s ∂b (i) h (t) = -2h (i) c (t) s∈S x(s)x T (s) x(s) 2 b (i) h . (9) Moving b (i)
h by a step length 1 2 λ(t) in the negative direction of this gradient, since λ(t) is small, we have:

b (i) h = I + λ(t)h (i) c (t) s∈S x(s)x T (s) x(s) 2 b ′ (i) h (10)
≈ s∈S I + λ(t)h (i) c (t) x(s)x T (s) x(s) 2 b ′ (i) h . (11)
This approximation amounts to a successive rotation of the basis for each component vector in the episode.

For stability of the solution, Kohonen proposed to multiply the learning rate λ(t) by x(s) xL i (s) and obtained the following slightly modified rotation operator:

M (i) c (t) = s∈S I + λ(t)h (i) c (t) x(s)x T (s) xLi (s) x(s) (12)
= s∈S P (i) c (s, t) . ( 13
)
That is how the operator

P (i)
c (s, t) was developed in the classical ASSOM basis updating.

In fact, the operator

P (i)
c (s, t) approximately maximizes another objective function as we shall discover now:

E m = i h (i) c s∈S xLi (s) x(s) P (X)dX . ( 14
)
Proof: The sample function of E m is:

E ms (t) = i h (i) c (t) s∈S xLi (s) x(s) . ( 15
)
Taking the gradient of E ms (t), we have

∂E ms ∂b (i) h (t) = h (i) c (t) s∈S ∂ xL i (s) x(s) ∂b (i) h . (16)
Moreover,

∂ xL i (s) 2 x(s) 2 ∂b (i) h = 2 xLi (s) x(s) ∂ xL i (s) x(s) ∂b (i) h . (17)
So,

∂E ms ∂b (i) h (t) = h (i) c (t) s∈S 1 2
x(s) xLi (s)

∂ xL i (s) 2 x(s) 2 ∂b (i) h . (18)
We have

h (i) c (t) s∈S ∂ xL i (s) 2 x(s) 2 ∂b (i) h (19) = ∂ l h (l) c (t) s∈S 1 -E s (t) ∂b (i) h (20) = - ∂E s ∂b (i) h (t) = 2h (i) c (t) s∈S x(s)x T (s) x(s) 2 b (i) h . (21) So ∂E ms ∂b (i) h (t) = h (i) c (t) s∈S x(s) xLi (s) x(s)x T (s) x(s) 2 b (i) h (22) = h (i) c (t) s∈S x(s)x T (s) xLi (s) x(s) b (i) h . ( 23
)
Taking a step λ(t) in the direction of this gradient, we get the rotation matrix

B (i) c (t) = I + λ(t)h (i) c (t) s∈S x(s)x T (s) xLi (s) x(s) . ( 24
)
When λ(t) is small, it is equivalent to [START_REF] Mcglinchey | Fast formation of invariant feature maps[END_REF].

M (i) c (t) as in

D. Further Boosting: Batch-mode Basis Vector Updating

Basis vector updating can be further boosted by working in a batch mode. We can avoid computing the value of xLi (s) in ( 6) by using the value computed previously during module competition. However, this could not be done inside the framework of FL-ASSOM since the subspaces are continuously changing in receiving each component vector of the episode. To save computation of xLi (s) , the batch-mode rotation operator B

(i) h = b ′ (i) h + ∆b (i) h , (25)
where ∆b

(i) h = s∈S α (i) c,h (s, t)x(s) . ( 26
)
The increment of each basis vector is thus a linear combination of component vectors in the episode. The difference between the updating rules (3) and ( 25) is that the former updates the basis vectors for each component vector one by one while the latter updates the basis vectors in a batch mode for the whole episode. The scalar parameter α

(i)
c,h (s, t) has the same form as [START_REF] Hase | Speech signal processing using Adaptive Subspace SOM (ASSOM)[END_REF]:

α (i) c,h (s, t) = λ(t)h (i) c (t) x T (s)b ′ (i) h xLi (s) x(s) . ( 27
)
The meaning of this equation is, however, a little different from that of [START_REF] Hase | Speech signal processing using Adaptive Subspace SOM (ASSOM)[END_REF]. Here in [START_REF] Mikolajczyk | An affine invariant interest point detector[END_REF] the basis vector updating is performed only after the whole episode has been received. Therefore, xLi (s) and x T (s)b

′ (i) h
can reuse the results previously calculated during module competition. What we need to do is only store the calculated values in registers and fetch them when needed. The computational load of ( 27) is thus trivial. Furthermore, the dissipation as well as orthonormalization of basis vectors can be performed only once for each episode without loosing accuracy since the basis vectors are not updated during the episode. The computational load can thus be further reduced. This method will be referred to as BFL-ASSOM for batch-mode fast-learning ASSOM.

Let us estimate the computational load of the BFL-ASSOM averaged on each component vector of the episode as we did for the basic ASSOM and for the FL-ASSOM. As has been mentioned, the calculation of α (i) c,h (s, t) according to [START_REF] Mikolajczyk | An affine invariant interest point detector[END_REF] needs only trivial computation. The majority of computation is in [START_REF] Harris | A combined corner and edge detector[END_REF]. Averaged on each vector in the episode, the computational load required by basis vector updating with the BFL-ASSOM is about HN multiplications and HN additions. Furthermore, since the dissipation and orthonormalization of basis vectors can be performed only once for each episode, the whole learning time can be further reduced.

E. Experiments

We first demonstrate by experiments that the BFL-ASSOM can generate topologically ordered invariant-feature filters as the basic ASSOM. The results of FL-ASSOM will be shown as the ground truth since the FL-ASSOM is mathematically equivalent to the basic ASSOM. Kohonen has shown that the ASSOM can generate basis vectors similar to Gabor filters for episodes subject to translation [START_REF] Kohonen | The Adaptive-Subspace SOM (ASSOM) and its use for the implementation of invariant feature detection[END_REF], [START_REF] Kohonen | Self-organized formation of various invariant-feature filters in the Adaptive-Subspace SOM[END_REF]. We shall show that the BFL-ASSOM is able to generate similar filters.

The input episodes are constructed from a colored noise image, which is generated by filtering a white noise image with a second-order Butterworth filter. The cut-off frequency is set to 0.6 times the Nyquist frequency of the sampling lattice. Each episode is composed of 6 vectors, each of which is formed on a circular receptive field composed of 349 pixels. The vectors in the same episode have only random translation of no more than 5 pixels in both the horizontal and the vertical directions. The episodes are generated on random locations of the colored noise image. The mean value of components of each input vector is subtracted from each component of the vector. In order to symmetrize the filters with respect to the center of the receptive field, the input samples are weighted by a Gaussian function symmetrically placed at the center of the receptive field with a full width at half maximum (FWHM) that varies linearly with respect to the learning step t from 1 to 16 sampling lattice spacings. Each vector is normalized to a unit vector before being sent to the ASSOM.

The ASSOM array is composed of 9 × 10 modules aligned in a hexagonal lattice with two basis vectors at each module. The basis vectors of all the modules are initially randomized and orthonormalized. The radius of the circular neighborhood function h (i) c (t) decreases linearly from 6.73 (= 0.5 × (9 2 + 10 2 ) 1/2 ) to 0.9 ASSOM array spacings with t. The learningrate factor λ(t) = 0.1•T /(T +99t), where T is the total number of learning steps and set to 30, 000 for the current experiment.

As shown in Fig. 3(a), the translation-invariant filters generated by the BFL-ASSOM and those by the FL-ASSOM are similar. The difference is only the different organization of the filters due to random initialization of the two networks. For both networks, the formed filters are similar to Gabor filters of different frequencies and different orientations. Moreover, filters of similar frequencies and orientations are formed at nearby sites. For either network, filters corresponding to b 1 and those corresponding to b 2 have the same frequencies at the same locations but 90 degrees of phase difference, which confirms orthogonality of the corresponding basis vectors. episode X = {x(s), s ∈ S}, where x(s) mean-subtracted and normalized, the projection error e(X) = s∈S x(s)x(s) 2 , where x(s) the projection of x(s) on the subspace of the winning module. e in the figure is the average of e(X) over all the training episodes. Fig. 3(b) confirms that the difference between the learning curve of the FL-ASSOM and that of the BFL-ASSOM is practically negligible.

In the second experiment, we compare the computational loads of the basic ASSOM, the FL-ASSOM and the BFL-ASSOM with respect to the input dimension N and the subspace dimension H. We record the elapsed CPU seconds for each method. The number of iterations are fixed to 1, 000. Each episode is composed of 6 vectors, which are generated randomly according to a uniform probability distribution. The rectangular ASSOM array contains 10 × 10 modules.

The timing results obtained by using C++ implementations are summarized in Table I, where the means and standard deviations of the elapsed CPU times after 20 runs with different initializations are recorded. As was anticipated, the basis vector updating time of the basic ASSOM increases rapidly with the input dimension and is the bottleneck of the learning procedure. With the FL-ASSOM, the basis vector updating time is dramatically reduced, increasing slowly with the input dimension and with the subspace dimension. It is no longer a bottleneck of the learning procedure. However, learning time outside basis vector updating is not reduced. Now with the BFL-ASSOM, the basis vector updating time is further reduced. Moreover, learning time outside the basis vector updating is also reduced considerably compared to the other two networks. Thus, the whole learning time decreases from 1, 956 seconds with the basic ASSOM to 16.2 seconds with the BFL-ASSOM when H = 4 and N = 400.

The relationship between the basis vector updating time and the input dimension or the subspace dimension for the three implementations of ASSOM is visualized in Fig. 4. The basis vector updating time increases approximately linearly with respect to the input dimension for the FL-ASSOM and for the BFL-ASSOM, but apparently nonlinearly for the basic ASSOM. In all the cases, the updating time increases approximately linearly with respect to the subspace dimension.

III. SALIENCY-BASED INVARIANT FEATURE CONSTRUCTION FOR IMAGE CLASSIFICATION

In this section, we explore a sequence of two schemes where the ASSOM is applied to saliency-based invariant feature construction for image classification. The implemented ASSOM may be the FL-ASSOM or the BFL-ASSOM. We will compare their performance in Section IV.

A. Salient-Point Single ASSOM Scheme (SPSAS)

The first scheme to be explored is based on a single ASSOM as shown in Fig. 5. Salient points are first detected from the input image by using the wavelet-based detector proposed in [START_REF] Laurent | In-depth analysis and evaluation of saliency-based color image indexing methods using wavelet salient features[END_REF]. Working with wavelets is justified by the consideration of the HVS for which multi-resolution, orientation and frequency analysis is of prime importance. In order to extract the salient points, a wavelet transform is firstly performed on the grayscale image. The obtained wavelet coefficients are represented by a zerotree structure [START_REF] Shapiro | Embedded image coding using zerotrees of wavelet coefficients[END_REF]. This tree is then scanned at a first time from leaves to the root to compute the saliency value at each node. Afterwards, this tree is scanned for the second time from the root to leaves in order to determine the salient path from the root to raw salient points on the original image. By working with grayscale images, the points located on boundaries of highlights or shadows are apt to be detected as salient whereas they are only caused by illumination conditions. To remove such false salient points, a gradient image is calculated by using the color invariants proposed by Geusebroek et al. [START_REF] Geusebroek | Color invariance[END_REF]. A threshold is then set to select the most salient points.

The ASSOM shall be trained to generate the appropriate feature filters based on the local regions (patches) around these salient points, which are supposed to carry essential information for image description and consequently for classification. The ASSOM can work on episodes composed of several component vectors, but construction of episodes is not necessary in a general situation if we want the ASSOM to learn subspaces. Although it is possible to construct episodes by artificially introducing transformations such as translation, rotation or scaling as in [START_REF] De Ridder | The Adaptive Subspace Map for image description and image database retrieval[END_REF], this process could generate artificial variants which might not really exist in the test set. So we prefer to use raw image patches to train the ASSOM.

Let p k , k ∈ {1, 2, . . . , K} be K salient points detected from the image I, which amount to K patches x k . These patches are fed into a single ASSOM with |I| modules, which was previously trained in an unsupervised way on all patches extracted from all categories of images in the training set. For each patch x k , the module i (i ∈ I) generates an energy xkLi 2 , with xkLi being the orthogonal projection of x on the subspace L i of the module i. Energies generated by all the modules construct an activity map, which is a vector We remark that this activity map construction corresponds to a soft (or fuzzy) clustering process in the sense that, in every module of the ASSOM, the input patch has a membership defined by the energy. It has been shown that a fuzzy clustering process has advantages over a hard one [START_REF] Lu | A method of face recognition based on fuzzy c-means clustering and associated sub-NNs[END_REF]. Activity maps of all patches extracted from I are then accumulated to form a cumulative activity map y I characterizing the image I, which is defined by:

a k = xkL1 2 • • • xkLi 2 • • • xkL |I| 2 T . ( 28
) Input image Wavelet-Based Salient Point Dectector Patch Extraction ASSOM Support Vector Machine Label Salient points { } K k p k , , 1 , K ∈ { } K k k , , 1 , K ∈ x

Cumulative

y I = K k=1 a k . (29)
This feature vector is then classified by an SVM previously trained on feature vectors of the training images. The implementation of the SVM used in this paper is part of a publicly available machine learning tool collection WEKA [START_REF] Witten | Data Mining: Practical Machine Learning Tools and Techniques[END_REF]. Here we use a Gaussian kernel G(u, v) = exp(-α u-v 2 ), where α = 0.02 in our experiments.

B. Salient-Point Multi-ASSOM Scheme (SPMAS)

The SPSAS does not make use of the fact that we have the label information for images in the training set. It might lead to map modules that mix up different categories features and confuse the (SVM) classifier. A better strategy would be to use a specific ASSOM for each category. This idea was explored in [START_REF] Zhang | Handwritten digit recognition by Adaptive-Subspace Self-Organizing Map (ASSOM)[END_REF] for recognition of handwritten digits and produced promising results. But the size of images in their case is very small (25×20 pixels), permitting a direct ASSOM learning. In our case, the images have much larger sizes and cannot be directly dealt with under their framework.

The SPMAS replaces the single ASSOM in the SPSAS with an array of ASSOMs, each one being trained for one category of image patches. Separate ASSOMs for different categories permit the system to learn the individual feature sets more precisely than a single ASSOM for all categories. Let J be the set of image labels (categories), c j = K k=1 a (j) k , j ∈ J be the cumulative activity map generated by the j-th ASSOM. The new feature vector y I is a joint cumulative activity map formed by combining the |J| cumulative activity maps:

y I = c T 1 • • • c T j • • • c T |J| T , (30)
as illustrated in Fig. 6. The feature vector y I is sent to the SVM for classification, as in the SPSAS. A similar single-SOM versus multi-SOM framework was proposed in [START_REF] Tan | Recognizing partially occluded, expression variant faces from single training image per person with SOM and soft kNN ensemble[END_REF] for face recognition. Each face image is firstly partitioned into non-overlapping sub-blocks corresponding to local feature vectors. A single SOM is trained for all classes of feature vectors or a separate SOM is trained for each class. Face images are then represented by sets of SOM weight vectors. A soft k-nearest neighbor ensemble method is proposed to identify unlabeled images. The single-SOM scheme and the multi-SOM scheme showed similar performance in their experiments. A block-to-block comparison is used in their research for face identification, which assumes good calibration between training faces and testing faces. It is not suitable for general-purpose image classification problems, where good calibration is not guaranteed. Also, representation of the feature vectors by the weight vectors may be noisy, a smoothed representation such as the activity maps may be more appropriate.

ASSOM 1 { } K k k , , 1 , K ∈ x Image patches ASSOM j ASSOM | | Joint Cumulative Activity Map y y                       ∑ ∑ ∑ = = = K k k K k j k K k k 1 ) ( 1 ) ( 1 ) 1 ( a a a M M                       ∑ ∑ ∑ = = = K k k K k j k K k k 1 ) ( 1 ) ( 1 ) 1 (

IV. EXPERIMENTAL RESULTS

A. Multi-Category Classification

In the first experiment, we evaluate our system in terms of multi-category image classification on the SIMPLIcity database 1 , which is part of the well known Corel database and has been used to test the SIMPLIcity content based image retrieval system in [START_REF] Wang | SIMPLIcity: Semantics-sensitive integrated matching for picture libraries[END_REF]. The database consists of ten categories including African people and villages (Afr), beaches (Bea), buildings (Bui), buses (Bus), dinosaurs (Din), elephants (Ele), flowers (Flo), food (Foo), horses (Hor), mountains and glaciers (Mou), each containing 100 images of 384 × 256 pixels. Some representative examples from each category are presented in Fig. 7. The images in each category are divided into two equal parts: 50 for training and the other 50 for testing.

We focus on RGB color features in the experiment. The image patches are circles with 597 pixels (about a radius of 14), which amount to 1, 791 = 597 × 3 dimensional vectors. The mean value of components of each input vector is subtracted from each component of the vector. The training steps are empirically set to 80, 000 for the ASSOM in the SPSAS and 30, 000 for each ASSOM in the SPMAS, simply because there are more training patches available for the single ASSOM in the SPSAS than for each ASSOM in the SPMAS. We use a rectangular ASSOM lattice. The choice between rectangular or hexagonal lattice does not seem to be critical according to our experiments.

1) SPSAS Versus SPMAS:

We first compare the SPMAS to the SPSAS based on the FL-ASSOM. The experimental results are summarized in Table II. Both schemes suffered from the overfitting problem due to lack of data. The SPSAS shows less overfitting than the SPMAS when the size of the ASSOM is small. But overfitting of the SPSAS inflates more quickly than that of the SPMAS when the size of the ASSOM increases. In this sense, the SPSAS is more sensitive to the number of parameters than the SPMAS. The SPMAS shows better accuracies than the SPSAS on both the training set and the test set across different configurations, which confirms advantages of using separate ASSOMs to learn features of different categories. On the test set, the SPMAS has an improvement of accuracy of 18.2%-28% over the SPSAS. We will stick to the SPMAS in the following experiments.

2) BFL-ASSOM Versus FL-ASSOM: Table III summarizes the classification accuracies of the SPMAS with the BFL-ASSOM and with the FL-ASSOM of various sizes. The number of neurons in the first layer of each ASSOM module is fixed to 2. This table shows that the performance of the BFL-ASSOM is a little better than the FL-ASSOM. The reason could be that the BFL-ASSOM is a more accurate learning process deduced from the corresponding objective function than the FL-ASSOM. Thus the local features could be better structured with the BFL-ASSOM than with the FL-ASSOM. The improvement of the classification accuracy is 1% -3.8% on the test set if we replace the FL-ASSOM with the BFL-ASSOM. Taking the learning speed into account, the BFL-ASSOM seems more attractive than the FL-ASSOM. According to Table III the performance of the BFL-ASSOMbased SPMAS is nearly optimal when the ASSOMs are of the size 10 × 10. We will stick to the 10 × 10 BFL-ASSOM-based SPMAS in the following experiments. 3) SPMAS Feature Filter Visualization: Fig. 8 shows the feature filters generated from the 10 categories of training images by using the 10×10 BFL-ASSOM-based SPMAS. Two neurons are implemented in the first layer of each module. In this way, each ASSOM learned two 10 × 10 lattices of basis vectors. In order to show the correspondence between the components of the basis vectors and the R, G, B components of the input patches, each three subsequent components of the basis vectors are grouped back to form an "RGB" pixel. The components are normalized to the range [0, 255] with the value 128 corresponding to a component 0. "Pixels"of the basis vectors are organized to form the same shape as the input patches, i.e. a circle.

Each basis vector can be rotated to the opposite direction without altering the spanned subspace. Thus we could rotate some of the basis vectors in order to get maps where neighboring basis vectors appear more similarly. For example, the orangish basis vectors could be turned to bluish in the b 1 lattice of the elephant category although we did not do that. This would not change the performance of the ASSOM and the consequent classifier.

As we mentioned previously in Section II-B, the ASSOM subspace learning process is similar to a recursive PCA. In the basis vector orthonormalization process, the first basis vectors are only normalized, and the second ones are orthogonalized with respect to the first ones. Thus the first basis vectors are likely to capture the first principal components of the input feature subspaces while the second ones are likely to capture the second principal components. There are some observable characteristics of features corresponding to the various categories of images. For flowers, the first basis vectors show a distinct red tone because most of the flowers, at least in the SIMPLIcity database, have a red tone. For buildings, both basis vectors do not show distinct colors. This is also the case for the dinosaurs since the dinosaur pictures are artificial models without distinct colors. In Fig. 9, we show three exemplary images and their joint cumulative activity maps, where gray levels of the modules in the ASSOMs are proportional to the energies of the respective modules. For the first two images, the ASSOMs corresponding to the correct image categories are more activated than other ASSOMs. However, this is not always so obvious. For example, the third image strongly activated both the "Bui" ASSOM and the "Bus" ASSOM. In fact, many bus images in the SIMPLIcity database have buildings in the backgrounds. It is not surprising to find out that the "Bus" ASSOM has learned the concept of some building parts. That is why we choose to keep the activity maps of different ASSOMs and use the SVM to make the decision.

4) Effects of the Subspace Dimension:

The classification results with the subspace dimension H = 4, 6, 8, 10 are summarized in Table IV. From H = 2 to H = 10, the average accuracy on the test set is improved for 7.8%, which suggests that higher subspace dimensions capture more precisely the variances of input patterns. The price to pay is that a higher subspace dimension involves a heavier computational load. Even though the complexity of the system seems to increase with higher subspace dimensions, the "curse of dimension" does not seem to appear. The reason could be that the length of the feature vector y I is constant with respect to the subspace dimension. Also, the SVM is less prone to overfitting than some other methods since it can limit the complexity of the model by the number of support vectors [START_REF] Duda | Pattern Classification[END_REF]. Fig. 10 shows some examples of the results from the test set. We also performed a 5 times 5-fold cross validation with H = 10, which shows that the classification accuracy can reach a mean value of 85.5% with a standard deviation of 2.6%.

According to [START_REF] Wang | SIMPLIcity: Semantics-sensitive integrated matching for picture libraries[END_REF] and [START_REF] Deselaers | Features for image retrieval: a quantitative comparison[END_REF], the best classification accuracy ever met across various features is 84.1% on the SIMPLIcity data set, which is worse than what we obtained, i.e. 86% on the test set or a mean value 85.5% with the 5 times 5-fold cross validation. Considering difficulty of the classification problem, the performance of our system is promising. It should be emphasized that the good results are obtained solely with raw image patch learning by ASSOM without calculating any descriptors beforehand. It shows that the system can learn patterns directly from the input signal itself.

5) SPMAS Confusion Matrix:

The confusion matrix of the SPMAS is shown in Table V. Most of the test images are correctly classified. The best classified category is dinosaur, where all the test images are correctly recognized. The worst classified is the beach category, where only 36 (72% ) images are correctly recognized. 4 (8%) images in the beach category are classified as elephant. Examination of the beach category shows that most images in this category have yellow sands and blue sky as backgrounds. Blue sky can often be found in the elephant category while the color of sand is often similar to that of soil in the elephant category. 8 (16%) mountains and glaciers images are classified as beach. This is because mountains can often be found in the beach category and they often share the blue sky as the background. In fact, an object could appear in a range of contexts and a salient point could be located at the border of two different but adjacent objects. When this happens, the context plays an important role. So higher level features, such as global features or semanticsbased features, should be helpful in further improving the classification accuracy. Also for some objects, such as sand in beach images and soil in elephant images, colors may not be enough for discriminating them whereas other descriptors such as texture should be more useful.

B. Adult Image Filtering

In the second experiment, we apply the SPMAS to adult image filtering. There are respectively 733 adult and 733 benign training images, 377 adult and 467 benign test images. Each BFL-ASSOM was trained with 200, 000 iterations. The subspace dimension is 2 and the topology of the ASSOMs is a 10×10 rectangular lattice. The image patches at salient points are circles of 597 RGB pixels. The trained BFL-ASSOMs are shown in Fig. 11, where b 1 of the adult images exhibits evident orientations, while b 2 embodies an obvious yellow tone. In Fig. 12, we show the exemplary joint cumulative activity maps of an adult image and a benign image. For the adult image, the "adult" ASSOM is obviously more activated than the "benign" ASSOM. However, for the benign image, the difference is not that obvious. This is related to the difficulty in learning "the rest of the world", which is practically impossible to be sufficiently sampled. Also, maybe we should wonder what the "benign" ASSOM has really learned from "the rest of the world". It has probably just learned a flat distribution, where nothing is specially interesting. We use the true positive (TP) rate and the false positive (FP) rate to describe the performance of adult image filtering. The TP rate is the proportion of correctly blocked adult images and the FP rate is the proportion of incorrectly blocked benign images. The SPMAS shows a TP rate of 89.1% with an FP rate of 13.9% on the training set. The F1 measure is 0.878 on the adult class and 0.874 on the benign class. On the test set, the TP rate is 90.2% and the FP rate is 13.1%. The F1 measure is 0.874 on the adult class and 0.892 on the benign class. Fig. 13 shows some classified examples from the benign test subset. The correctly classified images cover a wide range of scenes, including people or objects with skin-like colors. The incorrectly classified images include people with exposed skin or non-human objects with large areas of skin-like colors. To deal with such false alarms, higher-level analysis of the scenes might be necessary, such as detection of humans and context analysis. The receiver operating characteristics (ROC) curve of our system is shown in Fig. 14. The area under the curve (AUC) is a high value 0.958.

We compared the SPMAS with some other adult image filtering systems based on skin detection. Jones and Rehg [START_REF] Jones | Statistical color models with application to skin detection[END_REF] proposed one of the best skin detectors in the literature with the Compaq database. They built an adult image detection system based on their skin detector. Five features are calculated from the skin detection output and two additional features correspond to the dimensions of the input image. A neural network is implemented for classification. The adult image filtering system proposed in [START_REF] Zheng | Blocking adult images based on statistical skin detection[END_REF] is based on a MRF (Markov random fields) skin detector. Nine features are extracted from the skin detection output, including the shape information of skin regions. A multi-layer perceptron (MLP) is then implemented for classification. These systems, as well as the proposed SPMAS, were evaluated on different databases. Despite this, the comparison still makes sense since all the databases are comprehensive and contain thousands of images randomly downloaded from Internet, which are likely to follow similar distributions. Jones and Rehg's system has a TP rate of 88.9% under an FP rate of 13.1%. The system proposed in [START_REF] Zheng | Blocking adult images based on statistical skin detection[END_REF] has a TP rate of about 87.1% under an FP rate of 13.7%. Apparently the SPMAS is competitive to these skin-detection-based adult image filtering systems.

In fact for the skin-detection-based systems, although detected skin provides evidence of adult content, the false "skin" or missed real skin may confuse the adult image detector as well. The SPMAS, however, does not depend on such a skin preprocessor which might be a source of false detection itself. It learns features directly from raw image patches. In this way, the SPMAS also saves considerable manual labor of labeling skin pixels for training of the skin detector.

V. CONCLUSIONS AND PERSPECTIVES

The ASSOM is useful for dimension reduction, invariantfeature generation and visualization. Our study reveals that the increment of each basis vector in the ASSOM learning is a scaling of the component vectors of the input episode, which leads to a fast implementation of the ASSOM, i.e. the FL-ASSOM. With the FL-ASSOM, the computational load of basis updating is linear to the input dimension, which was quadratic with a naive implementation of the traditional ASSOM. We discovered the objective function approximately maximized by the traditional ASSOM and further proposed a batch-mode fast implementation of the ASSOM, i.e. the BFL-ASSOM, where the increment of each basis vector is a linear combination of the component vectors in the input episode. Computational load can be further saved with the BFL-ASSOM. The algorithms previously proposed in the literature could be faster than the traditional implementation of the ASSOM, but not than the fast implementations proposed here. The acceleration of ASSOM learning is especially meaningful for images, which are usually associated with high dimensions.

Experimental results revealed superiority of the SPMAS to the SPSAS in saliency-based invariant feature construction image classification. In the SPMAS, one ASSOM is trained for each image category to make use of the fact that image labels are known for the training set. The feature vector of the input image is built by combining cumulative activity maps of different ASSOMs. The SPMAS showed promising performance on a 10-category image classification problem and on adult image filtering. Compared to other skin-detectionbased adult image filtering systems, one important advantage of the SPMAS is that the manual labor related to preparation of the skin training set can be saved.

There could be several other ways to improve the ASSOM or the subsequent SPMAS. For example, a non-uniformly distributed input signal will lead to ASSOM modules with subspaces too close (when these modules are in the dense zones of the input space) or too far away from each other (when they are in the sparse zones of the input space). The similar problem is often encountered in the SOM. Hierarchical SOMs have been proposed to mitigate this problem. Such map structure permits multi-resolutional representation of the input signal space and fast locating of the winner as demonstrated by Liu et al. [START_REF] Liu | Motif discoveries in unaligned molecular sequences using self-organizing neural networks[END_REF]. A similar strategy could be used for the ASSOM. In this paper, the ASSOMs are only trained on the RGB color feature. We can as well train different ASSOMs on different features, e.g. color and texture, and then combine these ASSOMs. This could be realized in a boosting framework (AdaBoost for example).

Fig. 1 .

 1 Fig. 1. (a) An ASSOM lattice with a hexagonal neighborhood. Each empty circle in the lattice represents a neural module as shown in (b). The gray region represents the neighborhood area of the winning module indexed by c at some learning step. (b) Neural architecture of a module in the ASSOM.

Fig. 2 .

 2 Fig. 2. An alternative view of the basis vector updating rule of ASSOM.

  (i) c (t) in (24) will be useful. As was done for the FL-ASSOM, by distributing b ′ (i) h to terms in the operator B

  (i) c (t), the basis vector updating rule becomes: b

Fig. 3 (Fig. 3 .

 33 Fig. 3. (a) The Gabor-like filters generated by the BFL-ASSOM compared to those by the FL-ASSOM on episodes subject to translation. b 1 : First basis vectors. b 2 : Second basis vectors. (b) Changing of the projection error e with the learning step t for the FL-ASSOM and for the BFL-ASSOM.

Fig. 4 .

 4 Fig. 4. The basis vector updating time with respect to (a) the input dimension at the subspace dimension H = 2 and (b) the subspace dimension at the input dimension N = 200. For sake of clarity, the updating time of the FL-ASSOM and that of the BFL-ASSOM are magnified by a factor of 50.

1 aFig. 5 .

 15 Fig. 5. SPSAS image classification architecture.

Fig. 6 .

 6 Fig. 6. Construction of the joint cumulative activity map in the SPMAS.

Fig. 7 .

 7 Fig. 7. Example images from the SIMPLIcity database.

b 1 b 2 Fig. 8 .

 28 Fig. 8. Feature filters for the 10 categories of the SIMPLIcity training set. b 1 : The first basis vectors. b 2 : The second basis vectors.

  THE BOLDFACED FIGURES CORRESPOND TO CORRECTLY CLASSIFIED TEST IMAGES True Classes ↓ Afr Bea Bui Bus Din Ele Flo Foo Hor Mou Afr 43

Fig. 9 .Fig. 10 . 1 b 2 b 1 b 2 Fig. 11 .Fig. 12 .

 9101221112 Fig. 9. Three images and their joint cumulative activity maps. Afr Bea Bui Bus Din Ele Flo Foo Hor Mou

Fig. 13 .

 13 Fig. 13. Classification results of some benign test images. Over the line: Correctly classified examples. Under the line: Incorrectly classified examples.

Fig. 14 .

 14 Fig. 14. ROC curve of the SPMAS for adult image filtering.

TABLE I TIMING

 I RESULTS OF THE BASIC ASSOM, THE FL-ASSOM AND THE BFL-ASSOM WITH RESPECT TO THE INPUT DIMENSION N AND THE SUBSPACE DIMENSION H . LEFT SUB-TABLE: THE BASIS VECTOR UPDATING TIME (VU). µ REPRESENTS THE MEAN VALUE AFTER 20 RUNS. σ REPRESENTS THE CORRESPONDING SAMPLE STANDARD DEVIATION; RIGHT SUB-TABLE: THE WHOLE LEARNING TIME (WL), WHICH INCLUDES THE TIME FOR MODULE COMPETITION, BASIS VECTOR UPDATING, BASIS VECTOR DISSIPATION AND ORTHONORMALIZATION. ALL THE TIMES ARE GIVEN IN SECONDS

			H=2	H=3	H=4			H=2	H=3	H=4
			µ VU	σ VU	µ VU	σ VU	µ VU	σ VU			µ WL	σ WL	µ WL	σ WL	µ WL	σ WL
	ASSOM	N =100 N =200 N =300 N =400	78.0 303 682 1,331	0.54 2.1 3.8 11 1,621 96.8 377 846	0.73 3.2 7.9 1,003 120 449 13 1,904	4.7 2.7 7.6 13	ASSOM	N =100 N =200 N =300 N =400 1,354 84.0 314 699	0.57 2.2 3.8 11 1,658 107 395 872	0.75 3.2 8.0 1,042 134 475 13 1,956	4.8 3.0 7.8 13
		=100	0.942 0.083	1.26 0.097	1.60	0.10		N =100	6.92 0.057	10.9 0.086	15.6	0.13
	FL	N =200 N =300	1.46 2.10	0.13 0.12	2.01 2.80	0.16 0.20	2.48 3.60	0.10 0.18	FL	N =200 N =300	12.5 18.5	0.14 0.13	19.9 29.2	0.19 0.23	28.4 41.9	0.26 0.43
		N =400	2.80	0.26	3.73	0.28	4.66	0.20		N =400	24.3	0.27	38.7	0.32	55.1	0.37
		N =100	0.488 0.059 0.699 0.071 0.972 0.071		N =100	2.50 0.044	3.62 0.032	4.82 0.037
	BFL	N =200 N =300	0.569 0.056 0.817 0.093 0.745 0.085 1.05 0.091	1.12 1.46 0.095 0.14	BFL	N =200 N =300	4.35 0.057 6.34 0.064	6.28 0.065 9.19 0.071	8.42 0.076 12.3 0.085
		N =400	0.930	0.10	1.46	0.10	1.84	0.12		N =400	8.39 0.093	12.1	0.10	16.2 0.075
	basis vector updating time (s)	600 900 1200 1500 1800	ASSOM FL-ASSOM BFL-ASSOM
		300
		0
		100 200 300 400
			input dimension

TABLE II COMPARISON

 II OF CLASSFICATION ACCURACIES BETWEEN THE SPSAS AND THE SPMAS, AVERAGED ON THE 10 CATEGORIES

	ASSOM	SPSAS	SPMAS
	size	Training	Test	Training	Test
	5 × 5	45.8%	42%	87%	70%
	10 × 10	71.2%	58%	100%	76.2%
	15 × 15	88.4%	56.8%	100%	77%

TABLE III COMPARISON

 III OF CLASSIFICATION ACCURACIES BETWEEN THE BFL-ASSOM AND THE FL-ASSOM, AVERAGED ON THE 10 CATEGORIES

	ASSOM	FL-ASSOM	BFL-ASSOM
	size	Training	Test	Training	Test
	5 × 5	87%	70%	88.2%	73.8%
	10 × 10	100%	76.2%	100%	78.2%
	15 × 15	100%	77%	100%	78%

TABLE IV CLASSIFICATION

 IV ACCURACIES OF THE 10 × 10 BFL-ASSOM-BASED SPMAS WITH DIFFERENT SUBSPACE DIMENSIONS. THE RESULTS ARE AVERAGED ON THE 10 CATEGORIES

		H = 2 H = 4	H = 6 H = 8	H = 10
	Training	100%	100%	100%	100%	100%
	Test	78.2%	78.6%	82.4%	84.6%	86%

TABLE V CONFUSION

 V MATRIX OF THE SPMAS.

http://wang.ist.psu.edu/∼jwang/test1.tar

ACKNOWLEDGMENT

The authors are grateful to the anonymous reviewers for their valuable comments and suggestions, which have helped to improve the quality of this paper.