Advances in Hyperspectral Image Classification\\ {[Earth monitoring with statistical learning methods]}
%\input{abstract.tex} %
01_introduction.tex
%!TEX root = fullpaper_R1.tex
\section*{\small{INTRODUCTION}} \label{sec:intro}

The technological evolution of optical sensors over the last few decades has provided remote sensing analysts with rich spatial, spectral, and temporal information. In particular, the increase in spectral resolution of hyperspectral images and infrared sounders opens the doors to new application domains and poses new methodological challenges in data analysis. Hyperspectral images (HSI) allow to characterize the objects of interest (for example land-cover classes) with unprecedented accuracy, and to keep inventories up-to-date. Improvements in spectral resolution have called for advances in signal processing and exploitation algorithms. This paper focuses on the challenging problem of hyperspectral image classification, %\footnote{In the remote sensing community, the term `classification' is often preferred instead of the (perhaps more appropriate) term `segmentation', which is typically used in machine learning and computer vision.}
which has recently gained in popularity and  attracted the interest of other scientific disciplines such as machine learning, image processing and computer vision. In the remote sensing community, the term `classification' is used to denote the process that assigns single pixels to a set of classes, while the term `segmentation' is used for methods aggregating pixels into objects, then assigned to a class. %In the paper, we will use this terminology. 

% Necessary to motivate why are HSI interesting to SPM readers! Signal characteristics!
One should question, however, what makes hyperspectral images so distinctive. Statistically, hyperspectral images are not extremely different from natural grayscale and color photographic images~(see chapter 2 of~\cite{CampsValls11mc}). Grayscale images are spatially {\em smooth}: the joint probability density function (PDF) of the luminance samples is highly uniform, the covariance matrix is highly non-diagonal, the autocorrelation functions are broad and have generally a $1/f$ band-limited spectrum. In the case of color images, the correlation between the tristimulus values of the natural colors is typically high. While the three tristimulus channels are equally smooth in generic RGB representations, opponent representations imply an uneven distribution of bandwidth between channels. Despite all these commonalities, the analysis of hyperspectral images turns out to be more difficult, especially because of the high dimensionality of the pixels, the particular noise and uncertainty sources observed, the high spatial and spectral redundancy, and their potential non-linear nature. Such nonlinearities can be related to a plethora of factors, including the multi-scattering in the acquisition process, the heterogeneities at subpixel level, as well as the impact of atmospheric and geometric distortions. These characteristics of the imaging process lead to distinct nonlinear feature relations, i.e. pixels lie in high dimensional complex manifolds. The high spectral sampling of HSI (of the order of one band each 5-10 $nm$ in the electromagnetic spectrum) also leads to strong collinearity issues. Finally, the spatial variability of the spectral signature increases the internal class variability. All these factors, in conjunction to the few labeled examples typically available, make HSI image classification a very challenging problem. As a result, the accuracy obtained with standard parametric classifiers commonly used for multispectral image classification is typically compromised when applied to HSI~\cite{CampsValls05}. 

Many of these limitations have been recently addressed under the framework of statistical learning theory (SLT)~\cite{Vapnik98}. SLT is a general framework for learning functions from data, which reduces to finding a \blue{linear function defined in a high (eventually infinite) dimensional Hilbert feature space $f\in{\mathcal H}$ that learns the relation between observed input-output data pairs $(x_1,y_1),$ $\ldots,$ $(x_\ell,y_\ell)\in{\mathcal X}\times{\mathcal Y}$, and that {\em generalizes} well.} 
%The function $f:{\mathcal X}\to{\mathcal Y}$ should assign a class label $y\in{\mathcal Y}$ to each input $x\in{\mathcal X}$.} 
Generalization is the capability of a method to extrapolate to unseen situations, i.e. the function $f$ should accurately predict the label $y^\ast\in{\mathcal Y}$ for a new input example $x^\ast\in{\mathcal X}$. %, i.e. $y^\ast = f(x^\ast)$. 
Generalization has recurrently appeared in  statistics literature for decades under the names of bias-variance dilemma, capacity control, or complexity regularization trade-off. The underlying idea is to constrain too flexible functions in order to avoid overfitting the training data. 

The SLT framework formalizes this intuition~\cite{Vapnik98} and seeks for prediction functions that optimize a \blue{functional ${\mathcal L}_{reg}$ that takes into account {\em both} an empirical estimation of the training error (loss), ${\mathcal L}_{emp}$, 
and an estimate of the complexity of the model (or regularizer), $\Omega(f)$:
$$
\begin{array}{lll}
{\mathcal L}_{reg} &=& {\mathcal L}_{emp} + \lambda~\Omega(f) \\[2mm]
& = &  \sum_{i=1}^\ell V(x_i,y_i,f(x_i)) + \lambda~\Omega(f) 
\end{array}
$$
where $V$ is a loss function acting on the $\ell$ labeled samples, and $\lambda$ is a trade-off parameter between the cost and the regularization.}
Different losses and regularizers can be adopted for solving the problem, involving completely different families of models and solutions. 
%Actually, this minimization functional is shared by methods rooted in SLT, such as regularized networks and support vector machines (SVM). The loss function should be appropriate for the problem at hand. Actually, note that monotonicity of $\Omega$ does not ensure a unique minimizer of ${\mathcal L}_{reg}(x,y,f(x))$. Therefore, one has to prevent the possibility of several solutions, and indeed 
\blue{To ensure unique solutions, many SLT algorithms use strictly convex loss functions. The regularizer $\Omega(f)$ limits the capacity of the classifier to minimize ${\mathcal L}_{emp}$ and favors smooth functions. } 
%Regularization is also tied to encoding prior knowledge (and invariances to specific signal characteristics) in the classification function.

The hyperspectral image processing community has contributed to the design of specific loss functions and regularizers to take the most out of the acquired images. For example, regularization appears explicitly in many HSI classifiers when trying to impose the spatial homogeneity of images, when including the wealth of user's labeling in active learning, or when exploiting the information contained in the unlabeled pixels to better describe the image manifold in semisupervised learning. %Encoding regularization or prior knowledge are actually two sides of the same token. Note that 
Classifiers should also be robust to changes in the image representation: small perturbations of pixels and objects in the image manifold should not produce big differences in the classification. This is why the inclusion of proper image representations and invariances is also an active field.

In this paper, we review the recent advances in HSI classification under the SLT framework. Section II presents  advances in HSI classification using active, semisupervised and sparse learning approaches. Section III summarizes the field of spatial-spectral regularization both in terms of feature extraction and advanced classifiers. Section IV covers the  field of adaptation of both classifiers and feature representations, and reviews solutions to encode invariances in HSI classifiers. We conclude in Section V and outline future challenges.

02_standard.tex

\section{Standard Regularized image classification}\label{sec:robust}

%Hyperspectral images are rich and informative, but also increase significantly the dimensionality of data. 
Before HSI, most of the models used in remote sensing were parametric, such as Gaussian Maximum Likelihood (GML) or Linear Discriminant Analysis (LDA). These methods, based on the estimate of the covariance matrix, were successful when dealing with early multispectral images, whose dimensionality was usually comprised between four and ten bands. HSI changed the rules, as the dimensionality of the pixel representation raised to hundreds. Standard parametric methods became either unfeasible or unreliable, since estimating the class-covariance matrices require many labeled samples, which are usually not available. For that reason, research turned to regularized approaches, trying to limit the model complexity in the ill-posed situation of high-dimension and limited number of samples.

\vspace{-0.25cm}\subsection{Regularized parametric approaches.} %A first way to proceed is to add a regularizer to one of the traditional parametric methods, to force the sample covariance matrix toward a simpler, more stable matrix. 
Hyperspectral covariance matrices are typically conditioned via Tikhonov's regularization, $\hat{\Sigma} = \Sigma + \lambda {\Omega}$. Standard ridge reduces to ${\Omega} = {I}$. If ${\Omega} = {\beta}^\top {I}{\beta}$, where ${\beta}$ are the coefficients of optimal scoring of samples, reduces to {\em penalized discriminant analysis}~\cite{Yu99}, which %\cite{Hastie95}
downweights samples whose regression coefficients vary locally. %Using ${\Omega} = {I} - {VV}^\top$, where ${V}$ is %an orthogonal %%% THIS MUST BE WRONG: IF V IS ORTHOGONAL VV'=I
%a projection matrix to a $k$ subspace, %reducing dimensionality from $d$ to $k$ dimensions, 
%corresponds to 
Alternatively, one can perform classification in a subspace of reduced dimensionality by using ${\Omega} = {I} - {VV}^\top$, where $V$ is a projection matrix to that subspace~\cite{Jensen10}. Other strategies based on committees of classifiers trained on spectral subsets of HSI have also been proposed to attain more robust covariance estimates in LDA models~\cite{Prasad08}. % and boost final performances through decision function of the spectral views on the problem~\cite{Prasad08}.  
Finally, standard and more advanced forms of regularization are reviewed in~\cite{Bandos09}.

\vspace{-0.25cm}\subsection{Regularized nonparametric approaches.} Although successful, parametric models still rely on the assumption of normality of the class conditional PDFs. Due to the complexity of HSI, such assumption rarely holds, which encouraged research towards non-parametric models, not making any assumption on the distribution underlying the data. %However, nonparametric models are related to a strong increase of the models' parameters (typically, there are at least as many parameters as there are datapoints) and therefore placed the search for regularized solutions as a priority. % THIS IS NOT AT ALL TRUE: ACTUALLY, ONLY KERNEL METHODS FOLLOW THIS STATEMENT; NNETS OR TREES OR BOOSTING OR TYPICALLY ARE SPARSE NONPARAMETRIC MODELS...
In the last decade, SVM has become the most widely used method in HSI classification research~\cite{CampsValls05,CampsValls09}. Unlike other nonparametric approaches, such as regularized RBF networks, SVM naturally implements regularization through the concept of maximum margin: maximizing the linear separability between classes is equivalent to minimize the $\ell_2$-norm of model weights %, $\w$, 
used as regularizer. %, $\Omega=\|\w\|_2^2$. 
%is sparse and regularized by nature, since it includes in its optimization a $\ell_2$ norm regularizer of the model parameters $\mathbf{w}$, minimizing the energy of the data used and forcing smooth solutions~\cite{CampsValls05}. The problem is composed by two competing terms
%$${\mathcal L}_{reg}(x,y,f(x)) = {\mathcal L}_{emp}(f(x),y) + \lambda~||\mathbf{w}||^2_l,$$
%\begin{align}
%\min_{\w,\xi} C\sum_{i=1}^n\xi_i +& ||\w||^2_2\\
%\nonumber \text{s.t. } & y_i(\langle\w,x_i\rangle+b)\geq 1-\xi_i \text{, } \forall i, \xi_i \geq 0
%\end{align}
%The left term corresponds to the minimization of the empirical error (loss) and the right term is the SLT regularizer controlling the flexibility of the model. Using Lagrange multipliers, the final decision function is of the form 
%\begin{equation}
%y^* = f(\x^*) = sign\Big(\sum_{i=1}^n\alpha_iy_i\langlex_i , \x^*\rangle + b\Big)
%\label{eq:predSVM}
%\end{equation}
%where $\alpha_i$ are the SVM coefficients and are non-zero only for a small amount of training points (called the \emph{support vectors}). %For this reason, SVM solution is sparse.

Despite the elegance and simplicity of a nonparametric solution, SVM still provides only a linear classifier. % of type $f(\x^*) = \langle\w,\x^*\rangle + b$. 
When confronted to HSI, the linear SVM is still suboptimal since %, even if regularized, 
it does not provide nonlinear boundaries among the classes. A solution to this problem was found by SLT: the structural risk minimization principle~\cite{Vapnik98} states that it is easier to learn with a simple classifier in a higher dimensional space $\mathcal{H}$, than working in the original space $\mathcal{X}$ with a complex method. The difficulty %to apply this principle 
resides in projecting the samples $x$ into the higher dimensional space $\mathcal{H}$ using a mapping function $\phi(\cdot)$ since  there are infinitely many possible functions. %To do this, a characteristic of the decision function of Eq.~\eqref{eq:predSVM} is used: the SVM solution only depends on pairwise dot products (similarities) between each training sample and the point to predict, $\langle x_i,x^*\rangle$. The feature vectors are not used as such. 
The reproducing kernel properties gives the solution here: one only has to express model weights as a linear combination of examples in ${\mathcal H}$, and then replace the appearance of the dot products between mapped samples with a reproducing kernel function $K(x_i,x^*) = \langle \phi(x_i), \phi(x^*)\rangle_{\mathcal{H}}$ that works with input data points and represents the similarity between the samples in $\mathcal{H}$. Hence one can solve nonlinear problems by solving optimal linear classifiers in an appropriate higher dimensional vector space {\em implicitly}. %the method becomes nonlinear and performs linear separation in $\mathcal{H}$, which corresponds to a nonlinear separation in $\mathcal{X}$. This is called the \emph{kernel trick} and made the success of SVM, 
The SVM then allows to perform nonlinear classification without increasing the number of parameters, and keeps the elegant regularized formulation of the linear case.

%As for the LDA case, taking advantage of the block-diagonal structure of the covariance matrix to construct an ensemble of weak SVM classifiers boosts the performance~\cite{Waske10}.




03_advanced.tex
%!TEX root = fullpaper_R1.tex
\section*{\small{ADVANCED REGULARIZED IMAGE CLASSIFICATION}}\label{sec:advanced}

%\vspace{-0.05cm}\paragraph*{Standard regularized image classification}\label{sec:robust}

Before HSI, most of the classifiers used in remote sensing were parametric, such as Gaussian maximum likelihood or linear discriminant analysis. These methods, based on the estimate of the covariance matrix, were successful when dealing with early multispectral images, whose dimensionality was usually comprised between four and ten bands. HSI changed the rules, as the increased dimensionality of pixels raised to hundreds. Standard parametric methods became either unfeasible or unreliable, since estimating the class-covariance matrices requires many labeled samples, which are usually not available. 
For that reason, research turned to include regularization, %, trying to limit the model complexity in the ill-posed problem of dealing with a limited number of high-dimensional samples. Approaches considered explicit regularization,  
either explicitly through Tikhonov's terms in the involved covariance matrices, or by perform classification in a subspace of reduced dimensionality~\cite{Bandos09,Jensen10}. 
%Hyperspectral covariance matrices are typically conditioned via Tikhonov's regularization, $\hat{\Sigma} = \Sigma + \lambda {\Omega}$. Standard ridge reduces to ${\Omega} = {I}$. If ${\Omega} = {\beta}^\top {I}{\beta}$, where ${\beta}$ are the coefficients of optimal scoring of samples, reduces to {\em penalized discriminant analysis}~\cite{Yu99}, which %\cite{Hastie95}
%downweights samples whose regression coefficients vary locally. %Using ${\Omega} = {I} - {VV}^\top$, where ${V}$ is %an orthogonal %%% THIS MUST BE WRONG: IF V IS ORTHOGONAL VV'=I
%a projection matrix to a $k$ subspace, %reducing dimensionality from $d$ to $k$ dimensions, 
%corresponds to 
%Alternatively, one can perform classification in a subspace of reduced dimensionality by using ${\Omega} = {I} - {VV}^\top$, where $V$ is a projection matrix to that subspace~\cite{Jensen10}. Other strategies based on committees of classifiers trained on spectral subsets of HSI have also been proposed to attain more robust covariance estimates in LDA models~\cite{Prasad08}. % and boost final performances through decision function of the spectral views on the problem~\cite{Prasad08}.  
%Finally, standard and more advanced forms of regularization are reviewed in~\cite{Bandos09}.

Although successful, parametric models make strong assumptions about the normality of the class conditional PDFs or about the linearity of the problem. Due to the complexity of HSI, such assumptions rarely hold, which encouraged research towards non-parametric and nonlinear models. %, not making any assumption on the distribution underlying the data. %However, nonparametric models are related to a strong increase of the models' parameters (typically, there are at least as many parameters as there are datapoints) and therefore placed the search for regularized solutions as a priority. % THIS IS NOT AT ALL TRUE: ACTUALLY, ONLY KERNEL METHODS FOLLOW THIS STATEMENT; NNETS OR TREES OR BOOSTING OR TYPICALLY ARE SPARSE NONPARAMETRIC MODELS...
\blue{Several approaches have been introduced in the last decade in the field of hyperspectral image classification: kernel methods and support vector machines (SVM)~\cite{CampsValls05}, sparse multinomial logistic regression~\cite{Li12}, {neural networks~\cite{Ratle10}} and Bayesian approaches like relevance vector machines~\cite{Mianji11} and Gaussian Processses classification~\cite{Bazi10}. Nevertheless, the SVM has undoubtedly become the most widely used method in HSI classification research~\cite{CampsValls09}.} Unlike other nonparametric approaches, such as regularized RBF neural networks, SVM naturally implements regularization through the concept of maximum margin: \blue{given a linear classification function $f(x) = w^\top x + b$, maximizing the linear separability between classes is equivalent to minimize the $\ell_2$-norm of model weights $w$ used as regularizer, $\Omega(f)=\|w\|_2^2$.} Nonlinearity is also implemented via reproducing kernels, which allows to work in high dimensional Hilbert spaces implicitly, while still resorting to linear algebra operations~\cite{CampsValls05}. 

\iffalse
Despite the elegance and simplicity of a nonparametric solution, SVM still provides only a linear classifier. % of type $f(\x^*) = \langle\w,\x^*\rangle + b$. 
%When confronted to HSI, the linear SVM is still suboptimal since %, even if regularized, it does not provide nonlinear boundaries among the classes. 
A solution to this problem was found by SLT: the structural risk minimization principle~\cite{Vapnik98} states that it is easier to learn with a simple classifier in a higher dimensional space $\mathcal{H}$, than working in the original space $\mathcal{X}$ with a complex method. The difficulty %to apply this principle 
resides in projecting the samples $x$ into the higher dimensional space $\mathcal{H}$ using a mapping function $\phi(\cdot)$ since  there are infinitely many possible functions. %To do this, a characteristic of the decision function of Eq.~\eqref{eq:predSVM} is used: the SVM solution only depends on pairwise dot products (similarities) between each training sample and the point to predict, $\langle x_i,x^*\rangle$. The feature vectors are not used as such. 
The reproducing kernel properties gives the solution here: one only has to express model weights as a linear combination of examples in ${\mathcal H}$, and then replace the appearance of the dot products between mapped samples with a reproducing kernel function $K(x_i,x^*) = \langle \phi(x_i), \phi(x^*)\rangle_{\mathcal{H}}$ that works with input data points and represents the similarity between the samples in $\mathcal{H}$. Hence one can solve nonlinear problems by solving optimal linear classifiers in an appropriate higher dimensional vector space {\em implicitly}. %the method becomes nonlinear and performs linear separation in $\mathcal{H}$, which corresponds to a nonlinear separation in $\mathcal{X}$. This is called the \emph{kernel trick} and made the success of SVM, 
The SVM then allows to perform nonlinear classification without increasing the number of parameters, and keeps the elegant regularized formulation of the linear case.
\fi

%%%%%%%%%%%%%%%%

The effectiveness of SVM rapidly turned out to be insufficient to exploit the rich information contained in HSI. In all methodologies that follow, the functional to be optimized will consider additional information %of the signal, denoted hereafter as \blue{$(x_1',\ldots,x_m')$}, 
such as the one contained in unlabeled samples, ancillary data, or distinct signal characteristics. Such heterogeneous information is typically included in the HSI classifiers through additional regularizers. % that information Hence the decision function for a , $f(x)$, will becoming  $f(x,x')$. 

\begin{figure}[b]
\begin{tabular}{ccc}
Supervised & After SSL & After AL \\
\includegraphics[height=2.2cm]{./Figures/MKL_house_SSL_1_crop}&
\includegraphics[height=2.2cm]{./Figures/MKL_house_SSL_4_crop}&
\includegraphics[height=2.2cm]{./Figures/MKL_house_SSL_AL_crop}\\
\end{tabular}
\caption{Regularization of models with unlabeled samples: (left) purely supervised solution, (middle) semisupervised solution exploiting low density areas, and (right) active learning solution, where three new samples are labeled by a user.}
\label{fig:SSAL}
\end{figure}

\vspace{-0.25cm}\paragraph*{\small Regularization with unlabeled samples.} Recently, researchers started to exploit the abundant unlabeled information contained in the image itself, and new forms of regularization and priors were introduced. This is the field of \emph{semisupervised learning} (SSL, \cite{Bruzzone06,CampsValls07,Tuia09b}, central panel of Fig.~\ref{fig:SSAL}), where the minimization functional is modified to take into account %geometrical neighbors of a pixel in the manifold.  
the structure of the hyperspectral image manifold. 
%The models presented so far consider only the labeled pairs $\{x_i,y_i\}_{i=1}^n$. However, a HSI also provides a large amount of unlabeled samples: $\{x_j\}_{j=1}^m$. $n+m$ is the totality of the pixels of the image. To take advantage of unlabeled samples (and of the manifold structure they represent), \emph{semi-supervised} (SS) 
Usually, semisupervised algorithms modify the decision function of the classifier by adding an extra regularization term $\Omega_u$ that acts on \blue{both labeled and unlabeled examples: 
$$
\begin{array}{ll}
{\mathcal L}_{reg} =& {\mathcal L}_{emp} + \lambda~\Omega(f) + \lambda_u~\Omega_u(f).
\end{array}
$$
Several strategies to design the regularizer have been presented. One may use the graph Laplacian %${\bf L} = \D- {\bf W}$
 as a metric on the predictions to build $\Omega_u$. % (${\bf W}$ is the proximity graph and $\D$ is a diagonal matrix with the row sums of ${\bf W}$). 
Since regularization is performed on a proximity graph, the assumption enforced is that decisions on neighboring pixels in the data manifold should be similar~\cite{Gomez08}. Another possibility considers  
%Taking SVM as an example, if $\Omega_u = \frac{1}{u^2}\mathbf{f}L\mathbf{f}$, where $L$ is the graph Laplacian and \textbf{f} are the predictions on unlabeled data, regularization is performed on a proximity graph, thus enforcing that decision function of neighbors on the data manifold is similar~\cite{Gomez08}; if 
%$\Omega_u=\sum_{j=1}^u \mathcal{L}_{emp}(f(x'_j))$, where $\mathcal{L}_{emp}(f(x'_j))$ is the SVM hinge loss. Here %= \max(0,1-|f(x_j^*)|)$ 
%solutions with many unlabeled samples within the SVM margin are penalized, 
regularizers that enforce wide and empty SVM margins~\cite{Bruzzone06}.} Other strategies deform the kernel function by changing the metric induced using the unlabeled samples%, i.e. $K(\x,\x^*) = \alpha K(\x,\x^*) + (1-\alpha)K^C(\x,\x^*)$, where $K^C$ is a kernel function depending on the cluster structure of data
~\cite{Tuia09b,GomezChova10}.

% WITH THE OMEGA COLUMN
\iffalse
\begin{table}[!t]
\begin{center}
\caption{Summary of semi-supervised algorithms ($x'_j$ : unlabeled samples, $\mathcal{L}_{emp}$: hinge loss, $F_i$: classification of sample $i$).}
\label{tab:ssl}
\small
\setlength{\tabcolsep}{2pt}
\vspace{-0.15cm}
\begin{tabular}{p{1.7cm}|p{1.5cm}|p{2.3cm}|p{2.6cm}}
\hline
Assumption&Model &  Idea&  $\Omega_u$ \\\hline\hline
Low-density &TSVM & Look for the emptiest margin&  $\sum_{j=1}^u \mathcal{L}_{emp}(f(x'_j))$\\
\hline
Manifold, \underline{no out of} \underline{sample}&Laplacian SVM & Pixels close in the input space are also close in the graph&  $\mathbf{f}^\top\mathbf{L}\mathbf{f}$ (nearby samples are mapped close together) \\
\cline{2-4}
&Label propagation &  Spread class information on the graph of labeled/unlabeled &  {\scriptsize $\frac{1}{2} \displaystyle{\sum_{i,j=1}^n} W_{i,j} \Big\Vert \frac{F_i}{\sqrt{D_{ii}}} -$} {\scriptsize $-\frac{F_j}{\sqrt{D_{jj}}}\Big\Vert^2$} (nearby samples are classified in the same class) \\
\hline
Cluster & Cluster kernel & Deform the similarity (the kernel), then run standard & Increase similarity if samples  in the same cluster  \\
\cline{2-2}\cline{4-4}
&Mean map kernel & SVM   & Increase similarity if samples are mapped  close to centers in $\phi$ \\\hline
\end{tabular}
\end{center}
\end{table}
\fi

% WITHOUT THE OMEGA COLUMN
\begin{table}[!t]
\begin{center}
\caption{Summary of semisupervised algorithms used in HSI classification.}
\label{tab:ssl}
\small
\setlength{\tabcolsep}{2pt}
%\vspace{-0.15cm}
\begin{tabular}{p{1.7cm}|p{1.5cm}|p{5cm}}
\hline
Assumption&Model &  Idea   \\\hline\hline
Low-density & TSVM~\cite{Bruzzone06} & Look for the emptiest margin \\
\hline
Manifold & Label propagation~\cite{CampsValls07} &  Spread class information on the graph of labeled/unlabeled (nearby samples are classified in the same class) \\
\cline{2-3}
& Laplacian SVM~\cite{Gomez08} & SVM hinge loss plus Laplacian eigenmaps for manifold regularization: Pixels close in the input space are also close in the graph (nearby samples are mapped close together) \\
\cline{2-3}
& SSNN~\cite{Ratle10} & Neural network trained with gradient descent replaces SVM, graph regularization with loss that forces similar pixels to be mapped closely and dissimilar ones to be separated\\
\cline{2-3}
\hline
Cluster & Cluster kernel~\cite{Tuia09b} & Increase the similarity measure (kernel) if samples fall in the same cluster, then run standard SVM \\
\cline{2-3}
& Mean map kernel~\cite{GomezChova10} & Increase similarity if samples are mapped close to centroids in Hilbert space, then run standard SVM \\\hline
\end{tabular}
\end{center}
\end{table}


We evaluate the performance of semisupervised algorithms in an AVIRIS image acquired over the Kennedy Space Center (KSC), Florida in 1996, with a total of $224$ bands of $10$~nm bandwidth with center wavelengths from 400-2500 nm. The data was acquired from an altitude of $20$ km and has a spatial resolution of $18$ m. After removing low SNR bands and water absorption bands, a total of $176$ bands remains for the analysis. The dataset originally contained $13$ classes representing the various land cover types of the environment. Many different marsh subclasses were merged in a `marsh' class, resulting into the following $10$ classes: `Water' ($761$ labeled pixels), `Mud flats' ($243$), `Marsh' ($898$), `Hardwood swamp' ($105$), `Dark/broadleaf' ($431$), `Slash pine' ($520$), `CP/Oak' ($404$), `CP Hammock' ($419$), `Willow' ($503$), and `Scrub' ($927$).  \blue{The high  dimensionality and number of classes and subclasses pose challenging problems for the classifiers, especially when few labeled examples are available.}

Figure~\ref{fig:ksc-classified} illustrates classification results for cluster kernels, probabilistic mean map kernel, %transductive SVM (TSVM), 
label propagation, Laplacian SVM (LapSVM), and semisupervised neural networks (SSNN). We used $\ell=200$ labeled pixels ($20$ per class) and $u=1000$ unlabeled pixels. LapSVM, cluster kernels and mean map kernels perform similarly, and all improve the results of the label propagation whose training was particularly difficult in this high dimensional setting. % for which PCA dimensionality reduction was necessary. 
More homogeneous areas and better classification maps are observed in general for the mean map and bag kernels, and particularly for the SSNN, which efficiently deals with complex marsh areas (south east area of the image) and cope with large scale datasets.

\begin{figure*}[t!]
\small
\begin{center}
\begin{tabular}{cccc}
{\bf RGB} & {\bf SVM} (81.11\%, 0.82) & {\bf Bag kernel} (83.44\%, 0.83) & {\bf Mean map} (85.21\%, 0.84)\\ 
%\hspace{-1.6cm}\includegraphics[width=5cm]{./Figures/ksc_rgb_tm_low.pdf} &
%\hspace{-1.6cm}
\includegraphics[width=3.5cm]{./Figures/ksc_original.pdf} &
\includegraphics[width=3.5cm]{./Figures/ksc_svm.pdf} & 
\includegraphics[width=3.5cm]{./Figures/ksc_bag.pdf} &
\includegraphics[width=3.5cm]{./Figures/ksc_meanmap.pdf}\\
%\hspace{-1.3cm} 
& {\bf Label prop.} (70.57\%, 0.64) & 
%{\bf TSVM} (79.66,0.81)  
{\bf LapSVM} (83.11\%, 0.83) & 
{\bf SSNN (87.89\%, 0.87)}\\ 
%\hspace{-1.6cm}
%\hspace{-1.3cm}
& \includegraphics[width=3.5cm]{./Figures/ksc_labelpropagation.pdf} & 
%\includegraphics[width=3.5cm]{./Figures/ksc_tsvm.pdf} & 
\includegraphics[width=3.5cm]{./Figures/ksc_lapsvm.pdf} &
\includegraphics[width=3.5cm]{./Figures/ksc_ssnn2.pdf}
\end{tabular}
\end{center}
\vspace{-0.5cm}\caption{RGB composition %along with the ground truth (GT) available, 
and  classification maps with SVM, bag kernels, probabilistic mean map kernel, %TSVM, 
label propagation, LapSVM, and SSNN for the KSC image ($\ell=200$, $u=1000$). Overall accuracy and kappa statistic are given in brackets.}
\label{fig:ksc-classified}
\end{figure*}

\begin{table}[!t]
\caption{{Summary of  active learning algorithms~\cite{Tuia11} ($c$: number of candidates, $p$: members of a committee of learners).}}%, $x'_i$: already selected candidates).}}
\label{tab:AL}
\small
\setlength{\tabcolsep}{1pt}
%\vspace{-0.25cm}
%\begin{tabular}{p{1.3cm}|c|c|c|c c|p{2.6cm}}
%\hline
% Heuristic&   Batch& Type&  Classifier& \multicolumn{2}{c|}{ Diversity}& Models to train\\
% & &&& $S$& $SV$s& per iteration\\\hline\hline
% EQB&  \checkmark& M& $All$& $\times$ & $\times$& $p$ models\\
% AMD&  \checkmark&  M&  $All$&  $\times$ &  $\times$&  $p$ models\\
%
% MS&    \checkmark&  B&   SVM &  $\times$ &  $\times$&  SVM\\
% MCLU&   \checkmark&  M&  SVM&  $\times$ &  $\times$&  SVM\\
% SSC&    \checkmark&  B&  SVM&  $\times$ &  $\times$&  2 SVMs\\
% cSV&   \checkmark&  B&  SVM&  $\thickapprox$ &  \checkmark& SVM + distances to SVs \\
% MOA&  \checkmark&  B&  SVM&  \checkmark &  $\times$&  SVM + distances to already selected pixels\\
% MCLU-ABD&  \checkmark&  M&  SVM&  \checkmark&  $\times$&  SVM + distances to already selected pixels\\
% MCLU-ECBD&   \checkmark&  M&  SVM&  \checkmark &  $\times$& SVM + nonlinear clustering of candidates\\
% hMCS-i&   \checkmark&  M&  SVM&  \checkmark &  \checkmark&  SVM + nonlinear clustering of candidates and SVs\\%$-$&  
%
% KL-max&   $\times$&  M&  $$ $p(y|\x)$&  $\times$&  $\times$&  $(c-1)$ models\\
% BT&    \checkmark&  M&   $ p(y|\x)$&  $\times$&  $\times$&  1 model\\\hline
%\end{tabular}
%
%\end{table}
\blue{
\begin{tabular}{p{1.3cm}|p{1.3cm}| p{1.9cm}|p{1.5cm}|p{2.2cm}}
\hline
 Criterion&     Classifier& Uncertainty& Diversity& Models to train\\\hline\hline
 EQB&  $All$& Agreement of a committee& $\hspace{0.75cm}\times$ &  $p$ models\\\hline
 AMD&   $All$& Agreement of a committee&  $\hspace{0.75cm}\times$ &  $p$ models\\\hline
 MS&      SVM & Distance to SVM margin& $\hspace{0.75cm}\times$ &  1 SVM\\\hline
% MCLU&     SVM& Distance to SVM margins&  $\hspace{0.75cm}\times$ & 1 SVM\\\hline
% SSC&     SVM& Classify support vectors&  $\hspace{0.75cm}\times$ &    2 SVMs\\\hline
 cSV&    SVM&  Distance to SVM margin&  Spectral distance to current SVs& 1 SVM + distances to SVs \\\hline
 MOA&    SVM& Distance to SVM margin & Angular differences   & 1 SVM + distances to already selected samples\\\hline
% MCLU-ABD& SVM& Distance to SVM margins &     \checkmark&    SVM + distances to already selected pixels\\
 MCLU-ECBD& SVM& Distance to SVM margin & Different cluster assignment & 1 SVM + nonlinear clustering of $c$ samples \\\hline
 %hMCS-i&     SVM&  Distance to SVM margins &  Different cluster assignment &  1 SVM + nonlinear clustering of $c$ samples\\\hline%$-$&  
 KL-max&    Prob. output& Divergence of PDF if adding the candidate & $\hspace{0.75cm}\times$&    $(c-1)$ models\\\hline
 BT&      Prob. output& Difference in posterior of most confident classes& $\hspace{0.75cm}\times$&    1 model\\\hline
\end{tabular}
}
\end{table}

\vspace{-0.25cm}\paragraph*{\small Regularization via user's interaction.} Another possibility to cope with small sample problems is to provide additional labeled examples. \blue{This is possible since HSI represent land surfaces, usually physically reachable or that can be displayed in an image processing software. Therefore, the new samples can be collected either by photointepretation of the images (only if the classes can be recognized on screen) or 
by organizing field campaigns.} However, since providing additional samples is costly, the samples to be labeled must be selected carefully. To this end, \emph{active learning} (AL~\cite{Tuia11,Cra12}, right panel of Fig.~\ref{fig:SSAL})  has gained popularity in the last years: rather than proceeding by random sampling or stratification (i.e. sampling according to a measure of the expected variability within a class), AL uses the outcome of the current model to rank the unlabeled pixels according to their expected importance for future labeling. The aim is to detect the most difficult (and diverse) pixels for the current classifier. The top ranked pixels are then screened by a human operator, who provides the labels, enlarging the training set. With the enlarged training set, a new improved classifier is built and the process is iterated. Since AL focuses on difficult areas, % (the others are already well handled by the current model), 
it boosts the performances with fewer samples than those required by random sampling. Figure~\ref{fig:AL} illustrates an example of active learning model regularization  for the specific task of adapting classifiers to multiple scenes. % (see Fig.~\ref{fig:AL}).
%figure~\ref{fig:SSAL}c illustrates the principle. 








\vspace{-0.25cm}\paragraph*{\small Regularization through sparsity promotion.} 
Despite their high dimensionality, hyperspectral pixels belonging to the same class typically 
lie in a low-dimensional subspace. This observation was exploited in local 
semisupervised classifiers, and has been recently used in sparse signal representations. 
Here, the assumption is that pixels can be represented accurately 
as a linear combination of a few training samples from a structured dictionary. 

\blue{Note the connection with SVMs that embed the dictionary
(the training samples) into a high dimensional feature space ${\mathcal H}$. The use of the hinge loss} in the 
SVM functional induces a sparse solution, i.e. few training examples %(called {\em atoms} in the 
%sparsity literature or {\em support vectors} in the kernel literature) 
are selected. \blue{Recently, sparse kernel methods have been presented, such as the 
kernel matching pursuit, the $\ell_1$-SVM, the kernel basis pursuit or the generalized LASSO. In all these, 
the dictionary functions are the kernels centered around the selected `support vectors'. 
Alternatively, in~\cite{Chen12}, several sparse kernel approaches have been presented with a 
different philosophy: the target pixel is the test pixel itself, not a similarity evaluation, and the 
dictionary is composed by the training pixels in the feature space.} \blue{{In this paper}, a Basis Projection (BP) approach is 
used to promote sparsity with $\Omega = \|{w}\|_1$, as a relaxation of the more computationally demanding problem induced by using the $\ell_0$-norm.
%Since imposing sparsity induces an NP-hard problem, 
To solve the BP problem, greedy algorithms, such as the orthogonal matching pursuit (OMP) and the subspace pursuit (SP) methods, can be used. }
The dictionary can be obtained off-line or from the same 
image. % by PCA-like approaches or spectral unmixing models. 
The classification can be additionally improved by incorporating the contextual information
from the neighboring pixels into the classifier (see the next Section). 

\begin{figure}[t!]
\begin{center}
\includegraphics[width=6cm]{./Figures/sparse_aviris.pdf}
\end{center}
\vspace{-0.5cm}\caption{Performance measure with the estimated Cohen's kappa statistic, $\kappa$, for different sparsity-promoting classifiers.}
\label{fig:sparse}
\end{figure}

% Experimental comparison: 
In the example illustrated in Fig.~\ref{fig:sparse}, we compare the performance of the linear (SP, OMP) and kernel (KSP, KOMP) sparse HSI classifiers introduced in~\cite{Chen12}. The sparseness factor was tuned for best performance. 
We also included a linear SVM and the $\nu$-SVM using an RBF kernel, in which the parameter $\nu \in (0,1)$ controls the degree of sparsity. 
The RBF $\sigma$ parameter was tuned by standard $10$-fold cross-validation. Figure~\ref{fig:sparse} shows the 
results for these methods and different number of training samples in the standard AVIRIS Indian Pines 
hyperspectral image ($220$ spectral channels and spatial resolution $20$ m, shown in Fig.~\ref{fig:seg}). 
\blue{This is the standard benchmark hyperspectral image which is used here to allow comparison with results in~\cite{Chen12}.} 
We split the data into a training set (20\% of the available labeled pixels) and a test set (80\%). We trained the classifiers for different rates \{1, 5, 10, 15, 20, 25, 30\}\% of the training set, and show results for the test set that remained constant. 
%The procedure was repeated $10$ times and average results are reported. 
Nonlinear methods show a much better performance 
over linear approaches. In the linear case, SP clearly outperforms the rest, but when the nonlinearity 
is included all methods perform very similarly. % and the $\nu$-SVM is the method that improves its performance the most.


\iffalse
Recently, including sparse-promoting priors in general-purpose classifiers has yielded important 
improvements. In the second case, inducing sparsity can be accomplished introducing the $\ell_0$-norm as 
regularizer (aka entropy measure). However, this norm imposes difficulties in 
optimization and one typically resorts to the $\ell_1$-norm as a surrogate solution. A 
plethora of measures for sparse classification is now available, combining different 
combinations of loss plus regularizer. For instance, the LASSO method combines the least 
squares error and the $\ell_1$-norm regularizer, and has been found to have connections 
non-negative {\em garotte}, wavelet soft thresholding in orthonormal settings, 
basis pursuit, survival models, logistic regression, and the recent field of 
compressed sensing.

\begin{table*}[t!]
\begin{center}
\begin{tabular}{ll|cc}
\hline
Method                  & Author, Year       & $\min_{\bf w}\{\cdot\}$   & s.t. constraints \\
\hline
Ridge regression (RR)   & Hastie, 1983  & $\|y - \sum_i^d w_i x_i\|^2$         & $\|{\bf w}\|_2 = \sum_i^d w_i^2\leq t$ \\
LASSO                   & Hastie, 1983  & $\|y - \sum_i^d w_i x_i\|^2$          & $\|{\bf w}\|_1 = \sum_i^d |w_i|\leq t$ \\
Non-neg Garotte (NNG)   & Breiman, 1993 & $\|y - \sum_i^d c_i\hat w_i x_i\|^2$  & $c_i\geq 0,~~~~\sum_i^d c_i\leq t$ \\
Soft Thresh. (ST)       & Donoho, 1994  & $\|y - \sum_i^d w_i x_i\|^2 + \lambda\min({\bf w},t)$ & -- \\
Basis Pursuit (BP)      & Chen, 2000    & $\|{\bf w}\|_1$                           & ${\bf X}^\ast{\bf w} = y$ \\
BP denoising (BPDN)     & Chen, 2001    & $\|y - {\bf X}^\ast{\bf w}\|^2 + \lambda \|{\bf w}\|_1$ & -- \\
Dantzig Selector (DS)   & Cand\`es, 2007& $\|{\bf w}\|_1$                      & $\|{\bf e}{\bf X}^\ast\|_\infty \leq \lambda$,~${\bf e}={\bf y}-{\bf X}^\ast{\bf w}$ \\
\hline
\end{tabular}
\end{center}
\end{table*}

\fi









04_spatial.tex
%!TEX root = fullpaper_R1.tex

\section*{\small{SPATIAL-SPECTRAL IMAGE CLASSIFICATION}}\label{sec:spat}


%%% OLD FIGURE (2cols)
\iffalse
\begin{figure}[t!]
\small
\begin{center}
\setlength{\tabcolsep}{4pt}
\begin{tabular}{ccc}
RGB Composition & Reference & SVM (81.01)\\
\includegraphics[width=2.6cm]{./Figures/ChIV-JAB/pavia_3color}&
\includegraphics[width=2.6cm]{./Figures/ChIV-JAB/Pavia_reference}&
\includegraphics[width=2.6cm]{./Figures/ChIV-JAB/Pavia_SVM}\\
EMP (80.71) & EMAP (89.89) & DBFE+EMAP (94.50)\\
\includegraphics[width=2.6cm]{./Figures/ChIV-JAB/Pavia_EMP}&
\includegraphics[width=2.6cm]{./Figures/ChIV-JAB/Pavia_EAPall}&
\includegraphics[width=2.6cm]{./Figures/ChIV-JAB/Pavia_EAPall_DBFE}\\
\end{tabular}
\caption{RGB composition along with the available reference data for the ROSIS-03 Pavia University Area data set (103 spectral channels and spatial resolution 1.3m). Classification maps are shown for SVM on the original image only against spatio-spectral classification with the EMP, EMAP and EMAP after feature extraction by DBFE (81 data channels). Overall accuracies in percentage for each method are reported in parentheses. 
%(In the text, the latex reference for this figure is $\backslash$ref\{fig:spat\})
}
\label{fig:spat}
\end{center}
\end{figure}
\fi

\begin{figure*}[t!]
\small
\begin{center}
\setlength{\tabcolsep}{2pt}
\begin{tabular}{ccccccc}
RGB Composition & Reference & SVM %&EMP  
& EMAP  & DBFE+EMAP  & EMAP+CSVM \cite{CampsValls06}  & EMAP+GCSVM \cite{Li13} \\
&& $\omega$ & $s$ & $s$ & $\omega + s$ & $\omega + s$\\
&& (81.01) %(80.71) 
& (89.89) & (94.50) & (97.80) & (98.09) \\
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/pavia_3color}&
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_reference}&
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_SVM}&
%\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_EMP}&
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_EAPall}&
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_EAPall_DBFE}&
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_gustavocrossOK}&
\includegraphics[width=2.3cm]{./Figures/ChIV-JAB/Pavia_gckcrossok}\\
\end{tabular}
\vspace{-0.5cm}
\caption{RGB composition along with the available reference data for the ROSIS-03 Pavia University area data set (103 spectral channels and spatial resolution 1.3m). Classification maps are shown for SVM on the original image only against spatio-spectral classification with the EMP, EMAP and EMAP after feature extraction by DBFE (81 data channels), composite kernels wth cross-kernels and SVM~\cite{CampsValls06}, and generalized composite kernels with multinomial logistic regression~\cite{Li13}. Overall accuracies [\%] are reported in parentheses ($\omega$: using spectral bands, $s$: using spatial filters from PCA). 
%(In the text, the latex reference for this figure is $\backslash$ref\{fig:spat\})
}
\label{fig:spat}
\end{center}
\end{figure*}

%Since they observe the surface of Earth, 
Hyperspectral images live in a {\em geographical manifold}, in the sense that spatially neighboring pixels carry correlated information and that images are usually smooth in the spatial domain~\cite{CampsValls11mc}[see Ch.~2]. Accounting for spatial smoothness (i) provides less salt-and-pepper  classification maps, (ii) reveals the size and shape of the structure the pixel belongs to, and (iii) allows to discriminate between structures made of the same materials, but belonging to different land-use types. Spatial regularization has been widely used to improve classification~\cite{Fauvel13}. %Two basic families of spatial-spectral classifiers are encountered in the literature: 1) feature extraction (or preprocessing) techniques, and 2) filtering (or postclassification) methods. We will focus on the first one.
The joint exploitation of both \emph{spectral} and {\em spatial} information considers that either the loss, the regularizer, or both depend on the spatial neighbourhood of a pixel~\cite{Plaza09}. 

%\vspace{-0.25cm}\paragraph*{\small Adding spatial regularizers to spectral classifiers}
%first are solutions integrating priors on spatial smoothness based on Markov Random fields (MRF). MRF modeling has been shown to perform well for the classification of remote sensing images. However, classical MRF modeling (\emph{e.g.}, Ising, Potts) suffers when used for high spatial resolution data: neighboring pixels are highly correlated and the standard neighbor system definition does not contain enough samples to be effective. Unfortunately, a larger neighbor system imposes intractable computational problems, thereby limiting the benefits of conventional MRF modeling. Furthermore, algorithms involving MRF-based strategies traditionally require an iterative optimization step. 
 
\vspace{-0.25cm}\paragraph*{\small Spatial feature extraction.} 
A simple yet effective way to regularize for spatial smoothness is to enrich the input space with features accounting for the neighborhood of the pixels. This is usually done by using moving windows or adaptive filters applied to the spectral bands. These filtered images are then used to learn the classifier. %Two families of filters are the most used in the. \\
%The first is texture. Using the same crisp neighbor set employed by MRFs, textural features (occurrence and co-occurrence) can be also extracted from images. Such features provide relevant information about the granularity of the surface.\\
%Many families of spatial filters exist and have been applied: 
Standard filters based on occurrence or co-occurrence, morphological operators, Gabor filters or wavelets decompositions generally provide significant improvements over purely spectral classifiers. Among them, morphological filters are the most promising. In~\cite{Benediktsson05}, filtering was performed at many scales and an \emph{Extended Morphological Profile} (EMP) was used for classification. Proceeding in a multiscale fashion enables the adaptive definition of the neighborhood of a pixel according to the structure it belongs to. Filtering in HSI is more challenging than in multispectral images, and one typically resorts to compute the EMP based on only a few Principal Components (PCs) using morphological reconstruction operators. %rather computed on a subset of the HSI bands obtained with a FE method. 
 All the features are then fed to a classifier, either alone~\cite{Benediktsson05} or combined with the original spectral information~\cite{Fauvel08}. Furthermore, feature selection~\cite{Tuia10mkl} or extraction~\cite{Benediktsson05} can be used to find the relevant features. Recently, connected tree-based morphological operators have been investigated for the analysis of HSI~\cite{Mur10}. %They are based on tree-based image representations and offer new possibilities for extracting spatial information. 
These so-called attribute filters extract thematic attributes of the connected components of an image which are thresholded according to their geometry (area, length, shape factors), or texture (range, entropy). The multiscale version, Extended Morphological Attribute Profile (EMAP), has been also introduced. %For classification of HSI data, the attribute filters have also been considered in a multiscale fashion, leading to the Extended Morphological Attribute Profile (EMAP). %Since there are a multitude of possibilities for such filters, the feature selection part becomes even more important than in multispectral image classification. \devis{}{[ADD SPECIFIC REFERENCE? Dalla mura et al, GRSL 2011]} % The EMAP can become very highdimensional. Therefore, feature extraction or feature selection can be of great importance prior to its use in classification. 
%Figure~\ref{fig:spat} illustrates the remarkable gain in performance when using feature extraction on the EMAP, in conjunction with a SVM classifier.




\iffalse
%Alternative nonlinear feature extraction methods have been also exploited Principal Components Analysis (PCA) remains the most widely used method for this purpose, research of effective nonlinear approaches, such as Independent Component Analysis, Kernel PCA,  Decision Boundary Feature Extraction (DBFE)  or subspace methods is a growing current.
%The EMP is based on extracting a series of morphological filters from the hyperspectral data that contain most of the spectral information by some dimension reduction method. The EMP was first proposed with Principal Component Analysis (PCA), but its computation has later been extended to other feature extraction techniques, such as Independent Component Analysis (ICA), Kernel PCA (KPCA), . 
 opposite morphological operators (ex: opening ($\gamma$) and closing ($\phi$)), an ensemble of $s$ scales and a series of $m$ principal components projections, the EMP is 
\begin{equation}
{\scriptsize 
EMP = \left[
\begin{array}{c}
MP_1\\
\vdots\\
MP_m
\end{array}
\right] 
=
\left[ 
\begin{array}{ccccccc}
\gamma_1^{s} & ... & \gamma_1^{1} & PC_1 &\phi_1^{1} & ... & \phi_1^{s} \\
\vdots &&& \vdots &&& \vdots\\
\gamma_m^{s} & ... & \gamma_m^{1} & PC_m &\phi_m^{1} & ... & \phi_m^{s} \\
\end{array}
\right] \nonumber}
\end{equation}
%Consider the $m$ first principal components extracted from the hyperspectral image with PCA. For all components, the MPs are built as series of morphological filters results with structuring element of increasing size. For example, when using opening (Op) and closing (Cl) filters, the MP on the first principal component (PC$_1$) extracted form image $I$ will be $\text{MP}_1(I) = [\text{Op}_1^{S_3}, \text{Op}_1^{S_2}, \text{Op}_1^{S_1}, \text{PC}_1(I), \text{Cl}_1^{S_1}, \text{Cl}_1^{S_2}, \text{Cl}_1^{S_3}]$, with $S_1 < S_2 < S_3$. Then they are stacked to construct the EMP: 
%$$\text{EMP}(I)=\big[\text{MP}_1(I),\ldots,\text{MP}_m(I)\big].$$ 

%Additionally, EMP may contains some of the original spectral information, selected by the FE algorithm, and some spatial information carried by the morphological operators. 
%The EMP was originally used as an input to a classifier and gave good classification accuracies. However, the EMP contains only a part of the spectral information from the data. To overcome this problem, data fusion has been considered by using both the EMP and the original hyperspectral image and combining them into a stacked vector. Furthermore, feature extraction can also be applied on both feature vectors and the extracted features are concatenated in one stacked vector and classified by a SVM classifier. It has been shown that SVM can suffer from the curse dimensionality if many features are irrelevant or redundant. However, the feature extraction can overcome the problem.
\fi

%\vspace{-0.25cm}\paragraph*{\small Markov Random fields (MRF).} 
In Fig.~\ref{fig:spat}, the approaches based on mathematical morphology (EMP and EMAP extracted from the first four principal components) % \devis{}{extracted from the first principal component}) 
are used for classification of ROSIS-03 data from an urban area in Pavia, Italy. \blue{We selected this image to illustrate the capabilities of several spatial-spectral classifiers since urban areas monitoring at VHR typically requires the extraction of directional, rotational and scale features from objects}. A significant improvement in terms of classification accuracies with respect to the spectral SVM was achieved by applying the EMAP with four different attributes on the image {(area and diagonal of the bounding box of connected components, moment of inertia and standard deviation, see~\cite{Mur10}).}  On the other hand, some redundancy was observed in the original 144--dimensional filter vector of EMAP. Therefore, Decision Boundary Feature Extraction (DBFE) was applied on it. After extraction with DBFE, the accuracies improved significantly ($+5\%$), thus confirming the importance of feature extraction routines. Finally, comparison with composite kernels (see discussion below in Sect.\ref{sect:adv}~\cite{CampsValls06,Li13}) yielded improved  classification accuracy, but with a strong change of response in the right part of the image, where much more soil is predicted. This is explained by the fact that these two approaches jointly use spectral and spatial information, and are thus closer to the original spectral SVM result (for which the right part was partially predicted as soil).


\vspace{-0.25cm}\paragraph*{\small Spatial-spectral segmentation.} 
Another approach for the inclusion of spatial information is through image \emph{segmentation}, %Segmentation methods partition the image into non-overlapping homogeneous regions according to a criterion of homogeneity (usually based on the intensity of the pixels or on the texture of the neighboring patch). Each region in the segmentation map becomes the spatial neighborhood of the pixels it includes. The size of the regions depends on the smoothness of the homogeneity criterion: large neighborhoods are selected for homogeneous regions, while smaller region include highly distinct pixels. 
typically using watershed, mean shift and hierarchical segmentation~\cite{Fauvel13}. 
After segmentation, a supervised scheme assigns the pixels in the segments to the classes. Two approaches are the mostly used: in the first, the regions are treated as input vectors in a supervised classifier. In the second regions are considered as basins to post-process the class memberships attributed by a pixel-based classifier within each segment. 

The reverse view on the problem is proposed in~\cite{Tarabalka10b}, where a supervised classifier is used to produce confidence values for each pixel. Then, pixels with maximal confidence are used as seeds for a region growing algorithm. In~\cite{Mun11} segmentation and classification are linked by user-provided labels: working with a hierarchical segmentation of the data, the labels provided are used to isolate coherent clusters, both spatially and thematically, thus ending with the good segmentation and the labels of the segments. The number of queries is minimized with active learning. % (see Section~\ref{sec:robust}).

\begin{figure}[t!]
\small
\begin{center}
\setlength{\tabcolsep}{4pt}
\begin{tabular}{cc}
%False color & Reference \\
RGB & SVM (78.2\%, 0.75)\\
\includegraphics[width=3.8cm]{./Figures/ChIV-JAB/IndianPine_3color}&
\includegraphics[width=3.8cm]{./Figures/ChIV-JAB/IndianPinesSVM}\\
%\includegraphics[width=3.8cm]{./Figures/ChIV-JAB/IndianPine_Reference}\\
Majority Voting (90.8\%, 0.90) & Markers (91.8\%, 0.91)\\
\includegraphics[width=3.8cm]{./Figures/ChIV-JAB/IndianPinesMV}&
\includegraphics[width=3.8cm]{./Figures/ChIV-JAB/IndianPinesMarker}\\
%\multicolumn{2}{c}{SVM}  \\
%\multicolumn{2}{c}{(78.2\%, 0.75)}  \\
%\multicolumn{2}{c}{\includegraphics[width=3.8cm]{./Figures/ChIV-JAB/IndianPinesSVM}}\\
\end{tabular}
\vspace{-0.3cm}
\caption{RGB composition of the standard %along with the available reference data for the 
AVIRIS Indian Pine data set (200 spectral channels and spatial resolution 20m). Classification maps are shown for spatio-spectral classificaiton with the segmentation and classification with majority voting and segmentation with markers against SVM on the original image only. Overall accuracies and the kappa statistic for each method are reported in parentheses~\cite{Fauvel13}. 
}
\label{fig:seg}
\end{center}
\end{figure}

In Fig.~\ref{fig:seg}, two approaches based on segmentation and classification with majority voting and markers are applied to a 220-bands AVIRIS dataset over Indian Pines (Indiana, USA). Significant improvement in terms of overall classification accuracies and kappa statistic were achieved over the pixel-based SVM classifier. Using markers provided the best accuracies: a $+1\%$ improvement over the simple majority voting and more than $+13\%$ over the traditional pixel-based SVM classifier. Furthermore, as seen in Fig.~\ref{fig:seg}, it is clear that the classification and segmentation approaches provide a significantly more uniform classification map when compared to the purely spectral SVM classification map.

\vspace{-0.25cm}\paragraph*{\small Advanced spatial-spectral classifiers.} \label{sect:adv}
The main problem with spatial-spectral feature extraction approaches is the possibly high-dimensionality of the feature vectors to feed the classifier. This was alleviated in~\cite{CampsValls06} where dedicated kernels for the spectral and spatial information were combined. The framework has been recently extended to deal with convex combinations of kernels through multiple-kernel learning~\cite{Tuia10mkl} and generalized composite kernels~\cite{Li13}. In both cases, however, the methodology still relies on performing an {\em ad hoc} spatial feature extraction before kernel computation. Other alternatives in the literature considered %exploration of the space of spatial filters to discover the relevant filters~\cite{Tui13}, 
the definition of graph kernels that capture multiscale higher-order relations in a neighborhood without computing them explicitly~\cite{CampsValls10b}, and the modification of the SVM to seek for the spatial filter that maximizes the margin~\cite{Flamary12}. 

%\iffalse
A final alternative is to include contextual information with Markov Random Fields (MRF), which naturally include a spatial term on class smoothness in the energy function. %For this reason, they are appealing for the classification of remote sensing data, and have been applied successfully to single band or multispectral imagery. 
However, in the high dimensional context of HSI, the standard application of the neighbor system definition makes the problem computationally intractable, % does not contain enough samples to be effective, while a larger neighbor system imposes intractable computational problems, thereby limiting the benefits of conventional MRF modeling, further exacerbated by the iterative optimization step usually necessary (for example, using simulated annealing).
and therefore recent works have focused on joining MRF spatial priors and discriminative models in HSI classification~\cite{Tarabalka10,Moser13}. An excellent review of MRF spatial-spectral methods can be found in~\cite{Mos13}.
%\fi


\begin{table}[!t]
\caption{Summary of  spatial-spectral algorithms}
\label{tab:spat-spec}
\small
\setlength{\tabcolsep}{2pt}

\begin{tabular}{p{1.7cm}|p{2.5cm}|p{3.89cm}}
\hline
Type of \white{AAA} Approach&Model &  Idea \\\hline\hline
Spatial filters extraction &
Co-occurrence & Extract texture based on statistics of pairs of pixels in a neighborhood\\
\cline{2-3}
&EMP& Multiscale mathematical morphology (based on size)\\
\cline{2-3}
&EMAP& Multiscale  mathematical morphology (variety of attribute types) \\
\hline
Spatial-spectral segmentation & Segmentation and classification based on majority voting &  All pixels are assigned to the most frequent class inside a segmented region \\
\cline{2-3}
&Segmentation and classification based on markers & Most reliably classified pixels are selected as ``region markers'' for segmentation\\
\cline{2-3}
&Semi-supervised hierarchical clustering tree & Returns both classification and confidence maps. Active learning used to select informative samples.\\
\hline
Advanced spatial-spectral  &Composite and multiple kernels & Balances between spatial and spectral information with dedicated kernels \\
%\cline{2-3}
%&Multiple kernel learning & Convex combination of dedicated kernels with different contextual information and filter size\\
\cline{2-3}
classification&Graph kernels &Takes into account higher order relations in each pixel neighborhood \\
\cline{2-3}
&MRF&  Markov Random Field Modeling (probabilistic)  \\

\hline
\end{tabular}

\end{table}

05_adaptation.tex
%!TEX root = fullpaper_R1.tex

\section*{\small ADAPTATION AND INVARIANCES}% classifiers and representations}

One of the greatest challenges of modern HSI classification is the adaptation of classifiers between acquisitions that differ either by the zone they represent and/or the acquisition conditions such as illumination, angle and season, among other effects. Adaptation is a central issue in HSI classification. For example, the increase in revisit time of recent satellites 	%temporal resolution of image acquisitions of modern sensors allows data users to acquire several scenes of the same area, thus 
has improved multitemporal analysis of scenes. 
%multiangular studies and accurate change detection 
Nevertheless, algorithms must be able to adapt to changing situations. Generally, the direct application of classifiers trained on one image to new images leads to poor results: even if the objects represented in the images are roughly the same, differences %in reflectance, illumination, and atmospheric conditions
in acquisition induce significant {\em local} changes in the probability distribution function (PDF). % instead of simple global linear changes. 
These changes must be modeled and introduced in the classifiers. The concept of {\em adaptation} can be implemented at the levels of image pre-processing, robust and invariant feature extraction, or in the design of the classification algorithm. %In the first case, the goal is to adapt the radiometric properties of the image to be classified to those of the portion used for the learning of the classification algorithm. In the second case, the goal is to extract a subset of features from the original spectral channels that are as invariant as possible, while maximizing the class discriminability. In the third case, the scope is to adapt the classifier parameters estimated on the training samples extracted from a portion of an image to the properties of the classes observed in a different portion of it. 
%These approaches can be used in a collaborative framework for optimizing the adaptation results. %In what follows we will review the literature and will give experimental evidence of the performance of advanced methodologies recently introduced.

\paragraph*{\small Preprocessing.} 

The pre-processing phase can address adaptation through the use of radiometric correction techniques applied to the images. {\em Absolute corrections} aim at transforming the radiance measured at the sensor into surface reflectance. {\em Relative calibration} techniques adapt the radiometric properties between portions of an image or between images. Generally, absolute correction techniques require additional ground reference data that in many cases are not available or are difficult to collect. Relative calibration methods are often considered as a pragmatic alternative for adaptation. Among these methods, we recall histogram matching, relative radiometric normalization of time series~\cite{Canty04}, and multivariate histogram matching~\cite{Inamdar08}. %, and linear regression based on pseudo-invariant objects (e.g. man-made structures).

Images can be casted as point clouds in a geometrical space endorsed with an appropriate distance measure. Such a view is quite convenient because it allows us to move from {\em image adaptation} to {\em manifold adaptation}. Recent methods have explicitly considered the distortions occurring between image manifolds. In~\cite{Pet11}, multitemporal sequences for each pixel were aligned based on a measure of similarity between sequences barycenters, thus consisting into a global measure of alignment. In \cite{Jun11}, spectra of the pixels are spatially detrended using Gaussian processes in order to avoid shifts related to geometrical differences or to localized class variability. A recent principled approach tries to match graphs representing the data manifolds~\cite{Tuia12}. There, the graphs of the two domains are matched using a procedure aiming at maximizing their similarity, while at the same time preserving the original structure of the graphs. 

\if 0
As an example, we illustrate the graph matching method in~\cite{Tuia12} for the correction of angular effects. We studied a set of two HyMap hyperspectral images (128 channels) taken over the Barrax test site in 1999 during the DAISEX99 ESA campaign. The flight was designed to enhance angular effects. Two overpasses were acquired over the same site at noon: one in North-South direction along the solar principal plane and the other in East-West direction along the orthogonal plane (Fig.~\ref{fig:graphmatching}). This flight configuration allows the observation of an interesting directional anisotropy effect called \emph{hot spot} in the east-west flight line at noon. %In~\cite{Bei01}, a method for compensating for the angular effects based on the modeling of the bidirectional reflectance distribution function (BRDF) was proposed. However, such a compensation also requires knowledge of the land covers at the surface, that were estimated using a classification result. 
Graph-matching only looks at images geometrically and hence it performs the correction \emph{without} prior knowledge of the surfaces on the ground and their BRDF characterization. 
The manifold shape of images with and without hot spot is similar. %, see Fig.~\ref{fig:daisexManifoldsComp}. 
However, hot pixels are affected by an increase of their reflectivity values, and interestingly they occupy a particular spectral subspace so spectral adaptation can be feasible. Figure~\ref{fig:graphmatching} shows the performance of the graph matching algorithm for hot spot correction. We show the difference maps between the original acquisition and those obtained by two different ways of defining the manifold structure: in both cases, the graph matching procedure strongly reduces the angular effects, reducing the RMSE between the target domain ${\bf X}^{\mathcal T}$ and the corrected image ${\bf X}^*$ from 5.98 (no adaptation) to 4.55 ($k$-means) and to 4.19 (fuzzy $c$-means (FCM)).

% SPM wants big nice pictures with long self-explanatory captions... ;)
\begin{figure}[t!]
\begin{center}
\scriptsize
\setlength{\tabcolsep}{1pt}
\begin{tabular}{ccc}
Flight pattern  & Hot spot image, ${\bf X}^{\mathcal S}$ & Destination image, ${\bf X}^{\mathcal T}$ \\
\includegraphics[width=2.8cm]{./Figures/barrax_cruce} &
\includegraphics[width=2.8cm]{./Figures/immOK} \hspace{-2cm}\includegraphics[width=1.6cm]{./Figures/manifold_t2_trocito.pdf} &
\includegraphics[width=2.8cm]{./Figures/immHS} \hspace{-2cm}\includegraphics[width=1.6cm]{./Figures/manifold_t1_trocito.pdf}\\
$\|{\bf X}^{\mathcal T} - {\bf X}^{\mathcal S}\| = 5.98$ & $\|{\bf X}^{\mathcal T} -{\bf X}^*_{k\text{-means}}\| = 4.55$ & $\|{\bf X}^{\mathcal T} -{\bf X}^*_{\text{FCM}}\| = 4.19$\\
\includegraphics[width=2.8cm]{./Figures/hotVstarget} &
\includegraphics[width=2.8cm]{./Figures/correctedVstarget_km_spat_half_100centr_20models} &
\includegraphics[width=2.8cm]{./Figures/correctedVstarget_fcm_both_angleVZ100_fifth} \\
\end{tabular}
\end{center}
\vspace{-0.25cm}
\caption{Illustration of graph matching algorithms in the DAISEX dataset hot spot correction problem~\cite{Tuia12}. The hot spot is characterized by a sharp increase in reflectance when the surface is viewed in the same direction as it is illuminated by the Sun. Top panel motivates the common problem in HSI analysis: image acquisition implies in practice different illumination conditions due to the flight pattern. Whenever acquisition is done perpendicular to the principal plane, a strong hot spot phenomenon occurs which maximizes angular effects, ${\bf X}^{\mathcal S}$. Whenever acquisition is parallel to the principal plane, the hot spot vanishes, ${\bf X}^{\mathcal T}$. We show also inset the data manifolds for the images with and without  hot spot. In light blue, pixels corresponding to the hot spot area, which are those that should be adapted between the domains. Adaptation is thus strictly necessary to transform a flight line acquisition with hot spot to follow the statistics and manifold structure of the hotspot-free image line. Bottom panel shows the differences between the hotspot-free image ${\bf X}^{\mathcal T}$ and either the original angular hot spot image ${\bf X}^{\mathcal S}$ or its automatic correction ${\bf X}^*$ with manifolds defined either with $k$-means or FCM.}
\label{fig:graphmatching}
\end{figure}

\fi


\paragraph*{\small Adapting the classifier.} 
\label{sec:adaptC}
Learning a transformation between domains may be insufficient to handle all the perturbing factors, so alternative approaches are concerned with the adaptation of the classifier itself. From a machine learning perspective, the problem of classifier adaptation is studied in the framework of transfer learning, and in particular of domain adaptation. Domain adaptation reduces to learning from data in a source domain ${\mathcal S}$ (e.g. a portion of an image) to extrapolate to a different target domain ${\mathcal T}$ (another portion of the image or to another image). The problem has been payed attention in HSI classification lately~\cite{Bruzzone10}. In this setting, source and target domains are assumed to share the same set of information classes (exceptions to this constraint in~\cite{Tui11d,Jun13}) and to follow similar (but not the same) class distributions. Domain adaptation problems in remote sensing have been mainly addressed with semisupervised techniques, which exploit the labeled samples from ${\mathcal S}$ and the unlabeled samples from ${\mathcal T}$ in order to derive a classification rule suitable for the target domain. The most recent developments in this sense consider %maximum likelihood classifiers~\cite{Bru01}, cascade classification algorithms~\cite{Bruzzone2002}, 
semisupervised and domain adaptation SVMs~\cite{Bruzzone10}, Gaussian processes~\cite{Jun11,Jun13} and the mean-map kernel methods~\cite{GomezChova10}. 

\begin{figure}[t!]
\centerline{\includegraphics[width=8.8cm]{./Figures/AL-adapt}}
\vspace{-0.25cm}
\caption{Use of active learning to adapt a maximum likelihood classifier~\cite{Tui11d}. Left: given an image and a set of reference pixels in a first source area (${\mathcal S}$), we want to classify another spatially disconnected target area (${\mathcal T}$), by adding labels chosen actively in the reference of ${\mathcal T}$. Center: learning curves for the active (solid red line) and random sampling (dashed green line), evolving between the extremes of a model without adaptation (black dot at 150 samples), and another actively adapted that uses 550 samples randomly selected from ${\mathcal T}$ (blue dashed line). Right: classification maps in the target domain ${\mathcal T}$. The color of the bounding box in each map refers to the legend in the central plot.}
\label{fig:AL}
\end{figure}

Recently, active learning has been also used for adaptation assuming that some samples (as few as possible) from the target domain can be labeled by the user and added to the existing training set (defined on ${\mathcal S}$) in order to adapt the classifier to the target domain ${\mathcal S}$~\cite{Tui11d,Per12}. %On the other hand, training samples in ${\mathcal S}$ that do not fit the target distribution are removed from the learning process. 
This makes the adaptation process more robust than in the case of semi-supervised learning at the cost of requiring  additional labeled samples. {Figure~\ref{fig:AL} illustrates this principle for a $102$-bands image of the city center of Pavia acquired by the ROSIS-03 airborne sensor. The whole image (only a portion is shown) has a size of $1400 \times 512$ pixels and spatial resolution is  $1.3$ m. Five classes of interest (buildings, roads, water, vegetation and shadows) are considered, and a total of $206,009$ labeled pixels are available. We explore the potential of migrating a classifier built on a source area $\mathcal{S}$ with as few labeled pixels as possible from the rest of the image $\mathcal{T}$. First, the model trained with $150$ labeled pixels randomly drawn from ${\mathcal S}$ was directly applied to ${\mathcal T}$ and yielded a $\kappa_{\mathcal S}$ of $0.67$. Noticeably, a classifier built on $550$ samples randomly selected on ${\mathcal T}$ reached $\kappa_{\mathcal S} = 0.84$. This suggests that different regions of the same scene follow very different statistics for classification. In order to improve the first classifier, we enlarged the first training set with labeled samples drawn from ${\mathcal T}$, either taken through random sampling (RS, green dashed line) or with active learning (AL, red solid line). Using AL allows to concentrate efforts in areas where the first model is suboptimal, so performance is improved with respect to RS. After $400$ queries (thus, a model using $550$ training samples in total), random sampling yields similar performance to a model using $550$ randomly drawn pixels ($\kappa_{\text{RS}}=0.84$), while active learning improves the results with $\kappa_{\text{AL}}=0.89$. \blue{To reach the performance of the model using 550 random pixels, AL requires only 120 active queries (thus a total of 270 samples in the model)}.}

\paragraph*{\small Extracting and encoding invariances in the classifier.} 


\blue{Image classifiers must be robust to changes in the data representation within each land cover class}. The property of such mathematical functions is called `invariance'. 
A classifier should be invariant to object rotations, to changes in illumination, to the presence of shadows, and to the spatial scale of the objects to be detected. 
Extracting robust features (invariants) for classification and domain adaptation has been traditionally pursued by looking at the spatial or the spectral signal characteristics. On the one hand, {\em scale invariants} aim to make classifiers invariant to perturbations of object scales. In HSI classification, a single spatial scale is typically suboptimal because different classes exhibit diverse sizes, shapes, and internal variations. Multiscale classification schemes may alleviate these problems. Also, following wavelet-based representations and exploiting SIFT descriptors, {\em translation invariants} have been recently explored. 
On the other hand, {\em spectral invariants} are considered the fundamental descriptors of object structure, and are commonly employed to characterize canopy structure. %~\cite{Huang07}. In~\cite{Marchant02},  
Spectral invariance to daylight illumination allows for example improving multitemporal image classification. % where relative radiometric normalization aims at obtaining invariant pixels. %Efficient linear transformations have been proposed in this field~\cite{Canty02}.

%Even if a proper feature extraction has been conducted, a classifier may suffer from the variability within each land cover class.
 %Therefore, the classifier should be invariant to rigid transformations such as object translations, rotations or scaling. 
Incorporating invariances in SVM can be achieved by designing particular kernel functions that encode local invariance under transformations, or to generate artificial examples for training to which the model must be invariant. \blue{In the following example, we consider the latter possibility, with the Virtual SVM (VSVM) method}, which has been successfully exploited to encode scale, rotation, translation and shadow invariance in HSI classification~\cite{Izquierdo13}. 

In Fig.~\ref{results2}, we illustrate the use of the VSVM encoding (spectral) shadow-invariance for image classification. We use the same data acquired by the ROSIS-03 optical sensor of the city center of Pavia (Italy) used in Fig.~\ref{fig:AL}. %The image has a size of $1400 \times 512$ pixels of $1.3$ m spatial resolution, and consists of $102$ spectral bands covering the range $0.43$ to $0.86 \mu$m. 
%A total of $5$ classes of interest (buildings, roads, water, vegetation and shadows) are considered, and $206,009$ pixels were labeled by photointerpretation. 
We perform patch-based classification using only $50$ training patches of size $w=5$. The classes to be detected are, as for Fig.~\ref{fig:AL}, buildings, roads, water, vegetation and shadows. 
Virtual support vectors (VSVs) were generated  according to the observed  exponential behavior of the ratio shadow/sunlit as a function of the wavelength~\cite{Izquierdo13}. %~\cite{Yamazaki09}. 
%We should note that only the SVs belonging to the class `shadow' were used to generate VSVs. 
\blue{Numerical results, as well as the zoom on a detail of the classification map, show that VSVM leads to more accurate results than the standard SVM:  encoding shadow invariance reduces misclassifications on the bridge area and an overall more homogeneous classification over flat areas (see for example the crossroads in the center of the image)}. %The obtained numerical and visual results confirm the benefits of the invariance encoding. %Note that the method admits encoding several invariances simultaneously.

\begin{figure}[!t]
\small
\begin{center}
\setlength{\tabcolsep}{1pt}
\begin{tabular}{ccc}
RGB & SVM (0.79$\pm$0.11) & {\bf VSVM (0.84$\pm$0.09)}\\
\includegraphics[width=2.8cm]{Figures/rgb.pdf}&
\includegraphics[width=2.8cm]{Figures/svm.pdf}&
\includegraphics[width=2.8cm]{Figures/isvm.pdf}\\
\end{tabular}
\end{center}
\vspace{-0.15cm}
\caption{Experiment of patch-based classification with the virtual SVM encoding shadow invariance. True-color composite (left), and classification maps using the standard SVM (middle) and the VSVM (right) obtained using $50$ training pixels. Results are shown in parentheses in the form of (mean$\pm$ standard deviation of $\kappa$ in 20 realizations).\label{results2}}
\end{figure}





06_conclusions.tex
%!TEX root = fullpaper_R1.tex

\section*{\small CONCLUSIONS AND DISCUSSION}

This paper reviewed and analyzed the recent developments in hyperspectral image classification. Even though HSI follow similar spatial, spectral and spatial-spectral image statistics to those conveyed by conventional photographic images, the hyperspectral signals impose additional challenges related to their high dimensionality and heterogeneity. Therefore, even though standard techniques in image processing and computer vision may be transported drectly, HSI impose important constraints to develop efficient and effective classifiers.  

The use of methods derived from statistical learning theory (SLT) has been a driving factor in recent years. SLT constitutes a proper framework to tackle the problems posed by hyperspectral remote sensing images, which typically involve scenarios with high dimensional data and few training samples. SLT permit to embed numerical regularization in nonlinear classifiers, and also to design alternative forms of conditioning and incorporation of prior knowledge.  Additionally, classification is often improved by including spatially-based and manifold-based regularizers. \blue{SLT also allows to design sparse methods able to work in relevant feature subspaces, where compact and computationally efficient method can be run}. Finally, the SLT framework has allowed to include prior knowledge in a very fruitful way: for example, classifiers can now incorporate spatial and spectral invariances that disentangle ambiguities present in land cover classification.

The field is moving fast, and attracts research from the computer vision and machine learning communities. \blue{New approaches are introduced regularly and permit to tackle new scenarios issued from high resolution imaging (e.g. multitemporal, multiangular), while learning the relevant features via robust classifiers. It should be also noted that, with upcoming satellites, efficient algorithms for dimensionality reduction before classification and fast/parallel computing solutions will be necessary to accelerate the interpretation and efficient exploitation of hyperspectral images.}



abstract.tex
%!TEX root = fullpaper.tex
Hyperspectral images show similar statistical properties to natural grayscale or color photographic images. However, the classification of hyperspectral images is more challenging because of the very high dimensionality of the pixels and the small number of labeled examples typically available for learning. \blue{These peculiarities lead to particular signal processing problems, mainly characterized by indetermination and complex manifolds. %For these reasons, recent hyperspectral image classification research has been strongly influenced by the field of statistical inference. 
The framework of statistical learning has gained popularity in the last decade. New methods have been presented to account for the spatial homogeneity of images, to include user's interaction via active learning, to take advantage of the manifold structure with semisupervised learning, to extract and encode invariances, or to adapt classifiers and image representations to unseen yet similar scenes. 
%The design of a suitable regularized functional has led to the successful solution of many hyperspectral image classification problems. 
This tutuorial reviews the main advances for hyperspectral remote sensing image classification through illustrative examples. }% for each particular problem setting.

fullpaper_Final.tex
% Template for ICASSP-2010 paper; to be used with:
%          spconf_modified.sty  - modified ICASSP/ICIP LaTeX style file, and
%          IEEEbib.bst - IEEE bibliography style file.
% --------------------------------------------------------------------------
\documentclass[10pt]{article}
\special{papersize=8.5in,11in}
\usepackage{spconf_modified,amsmath,epsfig}

% Example definitions.
% --------------------
\usepackage{epsfig,cite}
%\usepackage{citesort}
\usepackage{amsmath}
\usepackage{amssymb}
\usepackage{color}
\usepackage{url}
\usepackage{multirow}
\usepackage{graphicx}
%Needed Occasionally:
\usepackage{epstopdf}
\newcommand{\red}[1]{\textcolor[rgb]{0.9,0,0}{#1}}   % Color roig
\newcommand{\green}[1]{\textcolor[rgb]{0,0.5,0}{#1}}   % Color verd
\newcommand{\blue}[1]{\textcolor[rgb]{0,0,0}{#1}}   % Color blau
\newcommand{\gray}[1]{\textcolor[rgb]{0.5,0.5,0.5}{#1}}  % Color gris
\newcommand{\white}[1]{\textcolor[rgb]{1,1,1}{#1}}  % Color gris

\newcommand{\devis}[2]{\textcolor[rgb]{0.5,0,0.9}{#2}}   % Color roig
\newcommand{\w}{\mathbf{w}}
\newcommand{\x}{\mathbf{x}}
\newcommand{\W}{\mathbf{W}}
\newcommand{\D}{\mathbf{D}}

\renewcommand{\baselinestretch}{1}

% Title.
% ------
\title{Advances in Hyperspectral Image Classification\\
{\small [Earth monitoring with statistical learning methods]}}
\name{Gustavo Camps-Valls, Devis Tuia, Lorenzo Bruzzone and J{\'o}n Atli Benediktsson}
\address{}


\begin{document}
\maketitle

%\begin{abstract}
%\input{abstract.tex}
%\end{abstract}

\input{01_introduction.tex}
\input{03_advanced.tex}
\input{04_spatial.tex}
\input{05_adaptation.tex}
\input{06_conclusions.tex}

\section*{\small ACKNOWLEDGMENTS}

\blue{This work was partially supported by the Spanish Ministry of Economy and Competitiveness (MINECO) under project LIFE-VISION TIN2012-38102-C03-01, and by the Swiss National Research Foundation under grant PZ00P2-136827.}

\section*{\small AUTHORS}

\blue{
\noindent {\bf\em Gustavo Camps-Valls} (gustavo.camps@uv.es) (SM'07) received a Ph.D. degree in Physics (2002) from the Universitat de Val\`encia, Spain, where he is currently an Associate Professor in the Electrical Engineering Dept. and head of the Image and Signal Processing group. His research is tied to machine learning for signal and image processing, with special focus on remote sensing data analysis. He has co-authored many journal papers and edited several books on kernel methods. The ISI lists him as a highly cited researcher. He is member of the Machine Learning for Signal Processing TC of the IEEE SPS, and serves as Associate Editor for IEEE Transactions on Signal Processing, IEEE Signal Processing Letters, and IEEE Geoscience and Remote Sensing Letters. Visit http://isp.uv.es/ for more information.\\
\noindent {\bf\em Devis Tuia} (devis.tuia@epfl.ch) received a Ph.D in environmental sciences at the Univeristy of Lausanne (Switzerland) in 2009. He was  then a postdoc researcher at both the University of Val\`encia, Spain and the University of Colorado at Boulder. He is now senior research associate at the laboratory of Geographic Information Systems, at the Ecole Polytechnique F{\'e}d{\'e}rale de Lausanne (EPFL). His research interests include the development of algorithms for information extraction and classification of very high resolution remote sensing images using machine learning algorithms. Visit http://devis.tuia.googlepages.com/ for more information.\\
\noindent {\bf\em Lorenzo Bruzzone} (lorenzo.bruzzone@ing.unitn.it) received the Ph.D. degree in Telecommunications Engineering from the University of Genoa, Italy, in 1998. He is currently a Full Professor at the University of Trento, Italy where he is the founder and the head of the Remote Sensing Laboratory. His current research interests are in the areas of remote sensing, radar and SAR, signal processing, and pattern recognition. He is the author (or coauthor) of many highly cited scientific publications in referred international journals, conference proceedings, and book chapters. He is the co-founder of the IEEE International Workshop on the Analysis of Multi-Temporal Remote-Sensing Images (MultiTemp) series. Since 2003 he has been the Chair of the SPIE Conference on Image and Signal Processing for Remote Sensing. Since 2013 he has been the founder Editor-in-Chief of the IEEE Geoscience and Remote Sensing Magazine. He is a member of the Administrative Committee of the IEEE Geoscience and Remote Sensing Society (GRSS) and a Distinguished Speaker of the IEEE GRSS. He is Fellow of IEEE. Visit http://rslab.disi.unitn.it for more information.\\
\noindent {\bf\em J{\'o}n Atli Benediktsson} (benedikt@hi.is) received the Ph.D. degree in electrical engineering from Purdue University, West Lafayette, IN in 1990. He is currently Pro Rector for Academic Affairs and Professor of Electrical and Computer Engineering at the University of Iceland. His research interests are in remote sensing, biomedical analysis of signals, pattern recognition, image processing, and signal processing. He was a co-founder of Oxymap, Inc. Prof. Benediktsson was the 2011-2012 President of the IEEE Geoscience and and Remote Sensing Society. He was Editor-in-Chief of the IEEE Transactions on Geoscience and Remote Sensing (2003-2008). He is a Fellow of IEEE and a Fellow of SPIE. Visit www.hi.is/$\sim$benedikt for more information.
}


\section*{\small REFERENCES}

\footnotesize{
\input{fullpaper_Final.bbl}}
%\bibliographystyle{IEEEbib}
%\bibliography{bibSPM-HSI,references}}

\end{document}

fullpaper_Final.bbl
\begin{thebibliography}{10}

\bibitem{CampsValls11mc}
G.~Camps-Valls, D.~Tuia, L.~G{\'o}mez-Chova, S.~Jimenez, and J.~Malo,
\newblock {\em Remote Sensing Image Processing},
\newblock Synthesis Lectures on Image, Video, and Multimedia Processing. Morgan
  and Claypool, 2011.

\bibitem{CampsValls05}
G.~Camps-Valls and L.~Bruzzone,
\newblock ``Kernel-based methods for hyperspectral image classification,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 43, no. 6, pp.
  1351--1362, 2005.

\bibitem{Vapnik98}
V.~Vapnik,
\newblock {\em Statistical Learning Theory},
\newblock Wiley, NJ, USA, 1998.

\bibitem{Bandos09}
T.~V. Bandos, L.~Bruzzone, and G.~Camps-Valls,
\newblock ``Classification of hyperspectral images with regularized linear
  discriminant analysis,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 47, no. 3, pp.
  862--873, 2009.

\bibitem{Jensen10}
A.~C. Jensen, M.~Loog, and A.~H.~Schistad Solberg,
\newblock ``Using multiscale spectra in regularizing covariance matrices for
  hyperspectral image classification,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 48, no. 4, pp.
  1851--1859, 2010.

\bibitem{Li12}
Jun Li, J.M. Bioucas-Dias, and A.~Plaza,
\newblock ``Spectral-spatial hyperspectral image segmentation using subspace
  multinomial logistic regression and {M}arkov random fields,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 50, no. 3, pp.
  809--823, 2012.

\bibitem{Ratle10}
F.~Ratle, G.~Camps-Valls, and J.~Weston,
\newblock ``Semisupervised neural networks for efficient hyperspectral image
  classification,''
\newblock {\em IEEE Trans. Geosci. Rem. Sens.}, vol. 48, no. 5, pp. 2271--2282,
  2010.

\bibitem{Mianji11}
F.A. Mianji and Ye~Zhang,
\newblock ``Robust hyperspectral classification using relevance vector
  machine,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 49, no. 6, pp.
  2100--2112, 2011.

\bibitem{Bazi10}
Y.~Bazi and F.~Melgani,
\newblock ``Gaussian process approach to remote sensing image classification,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 48, no. 1, pp.
  186--197, 2010.

\bibitem{CampsValls09}
G.~Camps-Valls and L.~Bruzzone,
\newblock {\em Kernel Methods for Remote Sensing Data Analysis},
\newblock J. Wiley \& Sons, NJ, USA, 2009.

\bibitem{Bruzzone06}
L.~Bruzzone, M.~Chi, and M.~Marconcini,
\newblock ``A novel transductive {SVM} for semisupervised classification of
  remote sensing images,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 44, no. 11, pp.
  3363--3373, 2006.

\bibitem{CampsValls07}
G.~Camps-Valls, T.~Bandos, and D.~Zhou,
\newblock ``Semisupervised graph-based hyperspectral image classification,''
\newblock {\em IEEE Trans. Geosci. Rem. Sens.}, vol. 45, no. 10, pp.
  3044--3054, 2007.

\bibitem{Tuia09b}
D.~Tuia and G.~Camps-Valls,
\newblock ``Semi-supervised remote sensing image classification with cluster
  kernels,''
\newblock {\em {IEEE} Geosci. Remote Sens. Lett.}, vol. 6, no. 1, pp. 224--228,
  2009.

\bibitem{Gomez08}
L.~G{\'o}mez-Chova, G.~Camps-Valls, J.~Mu{\~n}oz-Mar{\'i}, and J.~Calpe,
\newblock ``Semi-supervised image classification with {L}aplacian support
  vector machines,''
\newblock {\em {IEEE} Geosci. Remote Sens. Lett.}, vol. 5, no. 4, pp. 336--340,
  2008.

\bibitem{GomezChova10}
L.~G\'omez-Chova, G.~Camps-Valls, L.~Lorenzo~Bruzzone, and J.~Calpe-Maravilla,
\newblock ``Mean map kernels for semi-supervised cloud classification,''
\newblock {\em IEEE Trans. Geosci. Remote Sens.}, vol. 48, no. 1, pp. 207--220,
  2010.

\bibitem{Tuia11}
D.~Tuia, M.~Volpi, L.~Copa, M.~Kanevski, and J.~Mu{\~n}oz-Mar{\'i},
\newblock ``A survey of active learning algorithms for supervised remote
  sensing image classification,''
\newblock {\em {IEEE J. Sel. Topics Signal Proc.}}, vol. 5, no. 3, pp.
  606--617, 2011.

\bibitem{Cra12}
M.~M. Crawford, D.~Tuia, and L.~H. Hyang,
\newblock ``Active learning: Any value for classification of remotely sensed
  data?,''
\newblock {\em Proc. IEEE}, vol. 101, no. 3, pp. 593--608, 2013.

\bibitem{Chen12}
Yi~Chen, N.M. Nasrabadi, and T.D. Tran,
\newblock ``Hyperspectral image classification via kernel sparse
  representation,''
\newblock {\em IEEE Trans. Geosci. Rem. Sens.}, vol. 51, no. 1, pp. 217--231,
  2013.

\bibitem{CampsValls06}
G.~Camps-Valls, L.~G{\'o}mez-Chova, J.~Mu{\~n}oz-Mar{\'{\i}},
  J.~Vila-Franc{\'e}s, and J.~Calpe-Maravilla,
\newblock ``Composite kernels for hyperspectral image classification,''
\newblock {\em {IEEE} Geosci. Remote Sens. Lett.}, vol. 3, no. 1, pp. 93--97,
  2006.

\bibitem{Li13}
J.~Li, P.~R. Marpu, A.~Plaza, J.~Bioucas-Dias, and J.~A. Benediktsson,
\newblock ``Generalized composite kernel framework for hyperspectral image
  classification,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, in press.

\bibitem{Fauvel13}
M.~Fauvel, Y.~Tarabalka, J.~A. Benediktsson, J.~Chanussot, and J.~C. Tilton,
\newblock ``Advances in spectral-spatial classification of hyperspectral
  images,''
\newblock {\em Proc. IEEE}, vol. 101, no. 3, pp. 652 -- 675, 2013.

\bibitem{Plaza09}
A.~Plaza, J.~A. Benediktsson, J.~Boardman, J.~Brazile, L.~Bruzzone,
  G.~Camps-Valls, J.~Chanussot, M.~Fauvel, P.~Gamba, A.~Gualtieri,
  M.~Marconcini, J.C. Tilton, and G.~Trianni,
\newblock ``Recent advances in techniques for hyperspectral image processing,''
\newblock {\em Remote Sens. Environ.}, vol. 113, no. Supplement 1, pp.
  S110--S122, 2009.

\bibitem{Benediktsson05}
J.~A. Benediktsson, J.~A. Palmason, and J.~R. Sveinsson,
\newblock ``Classification of hyperspectral data from urban areas based on
  extended morphological profiles,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 43, no. 3, pp.
  480--490, 2005.

\bibitem{Fauvel08}
M.~Fauvel, J.~A. Benediktsson, J.~Chanussot, and J.~R. Sveinsson,
\newblock ``Spectral and spatial classification of hyperspectral data using
  {SVM}s and morphological profiles,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 46, no. 11, pp. 3804
  -- 3814, 2008.

\bibitem{Tuia10mkl}
D.~Tuia, G.~Camps-Valls, G.~Matasci, and M.~Kanevski,
\newblock ``Learning relevant image features with multiple kernel
  classification,''
\newblock {\em {IEEE Trans. Geosci. Remote Sens.}}, vol. 48, no. 10, pp. 3780
  -- 3791, 2010.

\bibitem{Mur10}
M.~Dalla~Mura, J.~A. Benediktsson, B.~Waske, and L.; Bruzzone,
\newblock ``Morphological attribute profiles for the analysis of very high
  resolution images,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 48, no. 10, pp.
  3747--3762, 2010.

\bibitem{Tarabalka10b}
Y.~Tarabalka, J.~Chanussot, and J.~A. Benediktsson,
\newblock ``Segmentation and classification of hyperspectral images using
  minimum spanning forest grown from automatically selected markers,''
\newblock {\em {IEEE} Trans. Sys. Man Cybernetics B}, vol. 40, no. 5, pp. 1267
  --1279, 2010.

\bibitem{Mun11}
J.~Mu{\~n}oz-Mar{\'i}, D.~Tuia, and G.~Camps-Valls,
\newblock ``Semisupervised classification of remote sensing images with active
  queries,''
\newblock {\em {IEEE Trans. Geosci. Remote Sens.}}, vol. 50, no. 10, pp.
  3751--3763, 2012.

\bibitem{CampsValls10b}
G.~Camps-Valls, N.~Shervashidze, and K.M. Borgwardt,
\newblock ``Spatio-spectral remote sensing image classification with graph
  kernels,''
\newblock {\em {IEEE} Geosci. Remote Sens. Lett.}, vol. 7, no. 4, pp. 741--745,
  2010.

\bibitem{Flamary12}
R.~Flamary, D.~Tuia, B.~Labbe, G.~Camps-Valls, and A.~Rakotomamonjy,
\newblock ``Large margin filtering,''
\newblock {\em IEEE Trans. Signal Proc.}, vol. 60, no. 2, pp. 648--659, 2012.

\bibitem{Tarabalka10}
Y.~Tarabalka, M.~Fauvel, J.~Chanussot, and J.~A. Benediktsson,
\newblock ``{SVM}- and {MRF}-based method for accurate classification of
  hyperspectral images,''
\newblock {\em {IEEE} Geosci. Remote Sens. Lett.}, vol. 7, no. 4, pp. 736
  --740, 2010.

\bibitem{Moser13}
G.~Moser and S.~B. Serpico,
\newblock ``Combining support vector machines and {M}arkov random fields in an
  integrated framework for contextual image classification,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 51, no. 5, pp.
  2734--2752, 2013.

\bibitem{Mos13}
G.~Moser, S.~B. Serpico, and J.~A. Benediktsson,
\newblock ``Land-cover mapping by {M}arkov modeling of spatial-contextual
  information,''
\newblock {\em Proc. IEEE}, vol. 101, no. 3, pp. 631--651, 2013.

\bibitem{Canty04}
M.~J. Canty, A.~A. Nielsen, and M.~Schmidt,
\newblock ``Automatic radiometric normalization of multitemporal satellite
  imagery,''
\newblock {\em Remote Sens. Environ.}, vol. 91, pp. 441--451, 2004.

\bibitem{Inamdar08}
S.~Inamdar, F.~Bovolo, L.~Bruzzone, and S.~Chaudhuri,
\newblock ``Multidimensional probability density function matching for
  pre-processing of multitemporal remote sensing images,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 46, no. 4, pp.
  1243--1252, 2008.

\bibitem{Pet11}
F.~Petitjean, J.~Inglada, and P.~Gancarski,
\newblock ``A global averaging method for dynamic time warping, with
  applications to clustering,''
\newblock {\em Pattern Recogn.}, vol. 44, no. 3, pp. 678--693, 2011.

\bibitem{Jun11}
G.~Jun and J.~Ghosh,
\newblock ``Spatially adaptive classification of land cover with remote sensing
  data,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 49, no. 7, pp.
  2662--2673, 2011.

\bibitem{Tuia12}
D.~Tuia, J.~Mu{\~n}oz-Mar{\'i}, L.~G{\'o}mez-Chova, and J.~Malo,
\newblock ``Graph matching for adaptation in remote sensing,''
\newblock {\em {IEEE Trans. Geosci. Remote Sens.}}, vol. 51, no. 1, pp.
  329--341, 2013.

\bibitem{Bruzzone10}
L.~Bruzzone and M.~Marconcini,
\newblock ``Domain adaptation problems: A {DASVM} classification technique and
  a circular validation strategy,''
\newblock {\em {IEEE Trans. Pattern Anal. Mach. Intell.}}, vol. 32, no. 5, pp.
  770--787, 2010.

\bibitem{Tui11d}
D.~Tuia, E.~Pasolli, and W.~J. Emery,
\newblock ``Using active learning to adapt remote sensing image classifiers,''
\newblock {\em Remote Sens. Environ.}, vol. 115, no. 9, pp. 2232--2242, 2011.

\bibitem{Jun13}
G.~Jun and J.~Ghosh,
\newblock ``Semisupervised learning of hyperspectral data with unknown
  land-cover classes,''
\newblock {\em {IEEE} Trans. Geosci. Remote Sens.}, vol. 51, no. 1, pp.
  273--282, 2013.

\bibitem{Per12}
C.~Persello and L.~Bruzzone,
\newblock ``Active learning for domain adaptation in the supervised
  classification of remote sensing images,''
\newblock {\em {IEEE Trans. Geosci. Remote Sens.}}, vol. 50, no. 11, pp.
  4468--4483, 2012.

\bibitem{Izquierdo13}
E.~Izquierdo-Verdiguier, V.~Laparra, L.~G\'omez-Chova, and G.~Camps-Valls,
\newblock ``Encoding invariances in remote sensing image classification with
  {SVM},''
\newblock {\em {IEEE} Geosci. Remote Sens. Lett.}, vol. 10, no. 5, pp.
  981--985, 2013.

\end{thebibliography}
