Classification of Multi-Frequency Polarimetric SAR Images Based on MultiLinear Subspace Learning of Tensor Objects

One key problem for the classification of multi-frequency polarimetric SAR images is to extract target features simultaneously in the aspects of frequency, polarization and spatial texture. This paper proposes a new classification method for multi-frequency polarimetric SAR data based on tensor representation and multi-linear subspace learning (MLS). Firstly, each cell of the SAR images is represented by a third-order tensor in the frequency, polarization and spatial domains, with each order of tensor corresponding to one domain. Then, two main MLS methods, i.e., multi-linear principal component analysis (MPCA) and multi-linear extension of linear discriminant analysis (MLDA), are used to learn the third-order tensors. MPCA is used to analyze the principal component of the tensors. MLDA is applied to improve the discrimination between different land covers. Finally, the lower dimension subtensor features extracted by the MPCA and MLDA algorithms are classified with a neural network (NN) classifier. The classification scheme is accessed using multi-band polarimetric SAR images (C-, Land P-band) acquired by the Airborne Synthetic Aperture Radar (AIRSAR) sensor of the Jet Propulsion Laboratory (JPL) over the Flevoland area. Experimental results demonstrate that the proposed method has good classification performance in comparison with the classic multi-band Wishart classifier. The overall classification accuracy is close to 99%, even when the number of training samples is small.


Introduction
Land cover classification is one of the most important applications of polarimetric synthetic aperture radar (SAR).The polarimetric SAR image can be classified more accurately than the single-channel SAR image, since polarization allows multi-channel measurements for the observed scene.At the same time, different frequencies are sensitive to different surface scales.The combined use of multi-frequency data can improve the classification accuracy, as illustrated by many researchers [1][2][3][4][5][6][7][8].However, it is not easy to automatically find effective target features for multi-band image classification.Some features may be effective at one frequency in distinguishing between targets, but they may fail at another frequency for these targets.Thus, it is crucial to find effective target features for classification of multi-frequency polarimetric SAR data.
In the past two decades, some classification schemes of multi-band polarimetric SAR data have been proposed.In view of the establishment of statistical characteristics of the scattering coherent matrix, Lee et al. [1] proposed a generalized maximum likelihood (ML) classifier, assuming that each frequency is statistically independent and that the probability density function (PDF) of the coherent matrix in single frequency follows the complex Wishart distribution.Assuming that the coherent matrix in two frequencies is statistically dependent, Famil et al. [2] introduced a entropy(H)/anisotropy(A)/alpha-Wishart classification scheme based on the complex Wishart density function of a 6 × 6 coherent matrix, which is constructed using the single-look complex data from the two frequency data.Datasets are classified with the improved Wishart classifier after an initial classification with the H/A/alpha classifier.Frery [3] presented several statistical distributions, including the Wishart distribution, K-distribution and G0-distribution, to model the statistics of different terrains of different frequencies, and the spatial context information was described with the multi-class Potts model.Then, the statistical information and contextual information are fused with the Bayesian framework.Gao et al. [4] put forward a multi-frequency classification algorithm by modeling the statistics of land with the mixture Gaussian or mixture Wishart PDFs, and then, an ML classifier is established for land covers.From the analysis of physical scattering characteristics of multi-band polarimetric SAR data, Freeman [5] proposed a hierarchical classification algorithm using different backscatter parameters based on scattering characteristics of different farmlands in different frequencies.Kouskoulas [6] estimated the PDF of principal scattering parameters of different terrains in different frequencies with the maximum entropy (ME) method and then classified short vegetation with the Bayesian hierarchical classifier.From the view of combining scattering features in different frequencies and then performing the classification with classic classifiers, Chen [7] presented a classification scheme using the dynamic learning neural network (NN) classifier to classify the combined features from multi-frequency.Lardeux [8] used the support vector machine (SVM) classifier to classify the high dimensional combined feature vector.
It has been demonstrated on the statistical analysis of multi-frequency polarimetric data that some scattering components are statistically correlated for single-frequency data, such as the HH and VV data being strongly correlated in homogeneous areas [9].The scattering characteristics are also correlated for different frequency data; they may be entirely relevant in some terrains, such as the bare soil at the L and C bands.In addition, objects usually have texture structures, and the scattering of one pixel may be correlated with its neighbors'.From the view of reducing the polarization, frequency and spatial correlation and extracting principal component of the data, Lee [10] developed a generalized principal component transform method to extract the principal component of multi-frequency polarimetric SAR images, and this method can maximize the signal-to-noise ratio and be tailored to the speckle noise characteristics.Trizna [11] presented a projection pursuit (PP) classification method, which projects scattering features of multiband polarimetric SAR images onto the subspace domain.Azimi-Sadjadi [12] extracted salient features of polarimetric data with principal component analysis (PCA) and then classified data with a neural network (NN) classifier.Chamundeeswari [13] used PCA to integrate the textural measures and wavelet components and classified the fusion features with the K-means classifier.Ainsworth [14,15] used two related nonlinear dimensionality reduction methods, i.e., locally linear embedding (LLE) and Isomap, to identify low dimensional structures of data; the subspace data are then classified with the geodesic distance.
Reducing the correlation of the multi-band data and extracting principal features either with PCA, PP, LLE or Isomap methods, the polarimetric and spatial features of all frequencies will be represented as a combined high dimensional feature vector.This will lead to the loss of the structural information of the original data, which can be separated into three different domains.In the last few years, a new principal feature extraction method was proposed based on the tools of multilinear algebra, by which not only the structure of the original data is maintained, but also effective principal features can be extracted.Many tensor decomposition methods summarized by Kolda [16] and low rank approximations of higher-order tensors, such as Lieven, which is presented in [17], have been applied to many fields, such as computer vision and signal processing.Using the tool of multilinear algebra, Lu et al. [18] extended PCA to the tensor data and proposed the multilinear principal component analysis (MPCA) method, such that the principal features of tensor structures can be extracted.Yan et al. [19] extended the linear discriminant analysis (LDA) to the tensor data and proposed the multilinear extension of linear discriminant analysis (MLDA) to extract features of tensor structures.Because the three-domain or even multi-domain features are simultaneously learned, the MPCA and MLDA methods are effective and useful when the features of the data are high dimensional and multi-domain.When these methods are used in target classification and recognition [18,20], the performance is better than that of the method of PCA or LDA.To possess the structure information and extract features more effectively, this paper proposes a new classification method for multi-frequency polarimetric SAR data based on tensor representation and multilinear subspace learning.MPCA and MLDA are combined to learn the multi-band polarimetric SAR data, which are represented by tensors.Firstly, each cell of the data is represented by a third-order tensor according to the frequency, polarization and spatial domains.Then, the MPCA algorithm is used to analyze the principal components of tensor data in all three domains, and MLDA is applied to improve the discrimination between different classes.Finally, the lower dimension subtensor feature is unfolded to a vector and then is classified by an NN classifier.
The rest of this paper is organized as follows.In Section 2, tensor algebra is summarized, and the MPCA and MLDA algorithms are described briefly.In Section 3, the proposed scheme is presented in detail.The classification results with multi-frequency polarimetric SAR data will be shown and analyzed in Section 4. A comparison will be also provided.Section 5 concludes the study.

Basic Tensor Algebra
The tensor algebra has been discussed by [16,17].A vector in a given basis is expressed as a one-dimensional array; a matrix is represented by a two-dimensional array; and a tensor with respect to a basis is represented by a multi-dimensional array.An n-th-order tensor is denoted as X = [x i 1 ...in ] ∈ K I 1 ×I 2 ×...×In , where n is the order of the tensor, which means that there are n path arrays.An element of X is denoted by x i 1 ...in , and I k is the dimension of the k-th path array.
The inner product of tensors A, B ∈ K A, A , and the corresponding distance of tensors A, B is A − B F .The n-mode product of tensor A ∈ K I 1 ×I 2 ×...×I N and matrix U ∈ K Jn×In are defined as is tensor unfolding in which the I n path index is kept constant, and the others paths are embedded into the I n × I n matrix sequentially.The n-mode product of tensor B = A× n U is equivalent to the n-mode unfolding of tensor Just as the singular value decomposition (SVD) of a matrix, any tensor has a similar higher order decomposition, named the higher order singular value decomposition (HOSVD) or Tucker decomposition.An I 1 × I 2 × ... × I n -tensor A can be expressed by an n-mode product, as follows.
where Σ is the core tensor of the decomposition.
• T denotes the matrix transpose operator.According to the equivalent equations of the n-mode product and the n-mode unfolding, Equation (1) can be represented as a Kronecker product of matrices by unfolding tensors A and Σ as follows.

MPCA
Similar to PCA for one-dimensional vector data, the MPCA of tensor data can be described as follows [18].
Supposing that {X m , m = 1, 2, ..., M } is a set of M tensor samples in K I 1 ×I 2 ×...×I N , the objective of MPCA is to search a set of projection matrices {U (n) ∈ K In×Jn , J n ≤ I n , n = 1, ..., N }, such that the energy of the new projected tensors is maximum, where U (n) are semi-orthogonal matrices, which project the original tensor space . The energy of new tensors can be denoted as , and the objective function is: Similar to the alternating iterative of HOSVD algorithm, Equation (3) can also be solved by an alternating iterative algorithm.For the given projection matrices Û(1) , ..., Û(n−1) , Û(n+1) , ..., Û(N) , supposing the tensor in W y is unfolding in the l-mode, then Û(l) is: where trace(.) is the trace of matrix and T .The quantity Equation ( 4) to be maximized can be recognized as a Rayleigh quotient problem, and columns of U (l) are the corresponding eigenvectors of the largest J l eigenvalues of matrix Ψ l .If the eigenvalues of Ψ l are denoted as λ 1 (t) , λ 2 (t) , ..., λ n (t) in the t-th iteration, is the minimum number when the energy ratio of output tensors to input tensors exceeds a threshold ρ, i.e., where ρ is usually set to 0.9∼0.999.The alternating iterative algorithm of MPCA can be described as below.Input elements: tensor samples X m ∈ K I 1 ×I 2 ×...×I N , m = 1, ..., M , the threshold ρ and the maximum number of iterations Γ.
Output elements: the projected matrices U (l) , l = 1, 2, ..., N and the output tensors Step 1: Calculate the mean of the tensor set Xm = 1 M M m=1 X m , and initialize N l=1 = I I l ×I l , where I I l ×I l is the identity matrix.
Step 2: In the t-th iteration, calculate l with the l-mode unfolding of tensor in W y , where Step 3: Calculate the eigenvalue decomposition l consists of the eigenvectors corresponding to the largest J (t) l eigenvalues using the criterion Equation (5).

MLDA
Similar to the LDA of one-dimensional vector data, the MLDA of tensor data can be described as follows [19,20].Supposing that {X m , m = 1, 2, ..., M } is a set of M tensor samples in K I 1 ×I 2 ×...×I N , the objective of MLDA is to search a set of projection matrices {U (n) ∈ K In×Jn , J n ≤ I n , n = 1, ..., N }, such that the inter-class distance is maximum and the intra-class variance is minimum, where U (n) are semi-orthogonal matrices, which project the original tensor space I 1 × I 2 × ... × I N into a tensor subspace J 1 × J 2 × ... × J N , and the projected tensors are . The objective can be given as follows, where C is the total class number, n c is the sample number in class c, Xc is the average tensor of the class c and X is the total average tensor of all samples.The numerator of Equation ( 7) is inter-class variances measured by the sum of the distances between Xc and X, and the denominator is intra-class variances measured by the sum of the distances between each tensor X i and its center tensor Xc i in class c i .Similar to the alternating iterative algorithm of MPCA, the optimization problem Equation ( 7) can also be solved by an alternating iterative algorithm.For given matrices Û(1) , ..., Û(n−1) , Û(n+1) , ..., Û(N) , unfolding the tensors in objective function with the l-mode, then: = arg max where T is the inter-class variance, T is the intra-class variance when the tensor is unfolded in l-mode and Φl has been defined in Equation ( 4).
The quantity Equation ( 8) is also a Rayleigh quotient problem, and the projected matrix U (l) consists of the eigenvectors corresponding to the most significant J l eigenvalues of the matrix S W (l) −1 S B (l) .
If the eigenvalues of S W (l) −1 S B (l) are λ 2 , ..., λ n at the t-th iteration, then J l can be chosen using criterion Equation (5).
The alternating iterative algorithm of MLDA is as follows.
Input elements: tensor samples X m ∈ K I 1 ×I 2 ×...×I N , m = 1, ..., M , the threshold ρ, the maximum number of iterations Γ and the class number C.
Output elements: the projected matrices U (l) , l = 1, 2, ..., N and the output tensors Step 1: Calculate the total average tensor Xm Step 2: In the t-th iteration, calculate S B (l) (t) and S W (l) (t) with the l-mode unfolding of the tensor, as follows. where Step 3: Calculate the eigenvalue decomposition S W (l) consists of the eigenvectors corresponding to the largest under the criterion constraint Equation (5).
Step 4: If the termination condition l ε, l = 1, ..., N or the number of iterations t = Γ is satisfied, then go to Step 5; otherwise, return to Step 2.

The Proposed Algorithm for Multi-Frequency Polarimetric SAR Data
Since the property of a cell in multi-frequency polarimetric SAR can be described in spatial, polarization and frequency domains, each cell of the image can be represented by a multi-order tensor.MPCA can be used to extract the principal component of the tensor data and meanwhile to possess the structure of the data in three different directions.This means that the principal component extraction can preserve the main characteristics of multiband data in three directions simultaneously.For the principal component of several kinds of targets, MLDA is useful to improve the discrimination between different classes.This means that we can also improve the discrimination in three directions for the multiband data.After extracting the important subtensor features.We can classify data by measuring subtensor features of different classes with the Frobenius distance of tensors.However, the subtensor of one class is of high-variability if the multiplicative speckle of data is strong.The distance measure is unstable for classification, which will result in poor classification.Therefore, we classify the subtensors with an NN classifier since the PDF of the subtensors is difficult to derive.The proposed algorithm is described as follows.

Tensor Representation of Multi-Frequency Polarimetric SAR Data
In the case of a single-frequency and single-look polarimetric SAR, each resolution cell is expressed by a 2 × 2 complex matrix named the Sinclair matrix, which contains scattering coefficients of a target in a pair of orthogonal polarization bases.

S =
S HH S HV S V H S V V (10) where H, V denotes the orthogonal polarization bases and the element S qp is the backscatter coefficient when the transmitting polarization is p and the receiving polarization is q.In a monostatic case, the reciprocity theorem holds, and the Sinclair matrix is symmetric, i.e.,S HV = S V H .Then, the matrix can be reduced to a three-dimensional scattering vector k.In the complex Pauli basis set, the vector becomes: For the multi-look case, the coherency matrix is often used, as follows.
where • H denotes the conjugate transpose; and L is the number of looks.
There is also a spatial structure for each class of targets; the scattering matrix of each pixel is correlated with its neighborhood pixels; one central cell can be represented by an N × N spatial cells around it.Thus, each cell of single-frequency polarimetric SAR data can be represented by a 3 × 3 × N × N tensor.If the number of frequency bands is M , then each cell can be represented by a fifth-order complex tensor C M ×3×3×N ×N .Since the N × N spatial cells are correlated with each other, we spread out the N × N spatial cell matrix to an N 2 -dimensional vector.Some elements of the T matrix are also correlated, especially the degrees of freedom being only five in the case of single-look data; thus, we spread out the 3 × 3 scattering matrix with a vector as (T 11 , T 22 , T 33 , Re(T 12 ), Re(T 13 ), Re(T 23 ), Im(T 12 ), Im(T 13 ), Im(T 23 )), where Re(•) is the real part of a complex element and Im(•) is its imaginary part.In this way, each cell of polarimetric SAR data in M bands can be represented by a third-order real tensor R 9×M ×N 2 .

Principal Tensor Component Extracting Based on MPCA and MLDA
After the multi-frequency data have been represented by third-order tensors, the total dimension of the data is large.The dimension of each cell is 9 × M × N 2 .Even if M and N are small, such as M = N = 3, the dimension is 243.If we spread out the tensor with a vector and then process it with a PCA or classify it with a classifier directly, the computation is huge.Since the principal component of the data is extracted based on the tensor in the MPCA method, the tensor is processed in three directions simultaneously.In this way, the dimension is small in each direction, and thus, the computation will be reduced greatly.In addition, if we process the data with vector representation, the structure information of the original data will be lost.Hence, we extract the principal subtensor component with the MPCA algorithm.After the principal component has been extracted, the MLDA algorithm is applied to improve the discrimination between different classes, and the dimension of the subtensor is reduced at the same time.The computation time of MLDA is also less than that of the LDA method, and the structure of the data in three different directions is still possessed.
The specific algorithm flow chart is shown in Figure 1.First, each cell of all of the multi-band data is represented by a third-order tensor R 9×M ×N 2 ; then, the training samples are processed with the MPCA algorithm, and the projection matrices and the projected tensors of training samples are outputted.After the subtensors are trained with the MLDA algorithm, we can get the final projection matrices and the final projected tensors of training samples.Finally, the principal components of the testing data are extracted with the projection matrices that have been trained.
Let U (l) N l=1 be the output projection matrices with the MPCA training algorithm, V (l) N l=1 be output projection matrices with the MLDA training algorithm and the testing tensors be X i , (i = 1, ..., M ), and then, the extracting subtensors are Z i , where: Figure 1.Algorithm flow chart of the subtensor extracting with the proposed algorithm.

Classification with Neural Network
After the subtensors have been extracted by the MPCA and MLDA algorithms, we can classify data using the Frobenius distance between testing subtensors and training subtensors.Since the measure of subtensors is unstable for the multiplicative speckle of the original data, the subtensor is spread out with a vector and classified with a neural network classifier.
Let the training sample be X i , i = 1, ..., M , the extracted subtensor be Z i , the corresponding vector be z i , the class labels be C i , C i ∈ {1, ..., N c }, the output of z i traversing the NN be y i and the corresponding desired output vectors be c i .Then, the training model of NN is: where e i = y i (w) − c i is the error of the i-th sample and w is the weight of the neural network.After the NN is trained and the weight of the neural network w is obtained, the output class label of the test sample X k is: where y k,j (w) is the value of the output node j for test sample X k .

Experimental Results and Analysis
Multi-frequency polarimetric SAR images acquired by AIRSAR with the P, L and C bands, over Flevoland in the Netherlands, are used to test the performance of the proposed classification algorithm.Pauli pseudo-color images of the three bands are presented in Figure 2a-c, and the ground-truth data are shown in Figure 2d, which is drawn according to the paper of Hoekman [21].There are 14 kinds of crops on this site, including potato, fruit, oats, and so on.Figure 2e is the types of crops and the corresponding colors in the ground-truth.
To demonstrate the effectiveness of the proposed method, three experiments are designed.In Experiment 1, we want to show the classification performance and results when the training samples are enough and the dimension of the subtensors is reasonable by setting high thresholds in both the MPCA and MLDA steps.The classification performance is compared with those of the classic Wishart classifier and the NN classifier for the tensor data, which are unfolded into a vector without dimension reducing.In Experiment 2, we want to compare the classification performance of the tensor learning method with the vector learning method, and we also want to compare the classification performances of different tensor learning methods, including the tensor learning only using MPCA and only using MLDA, as well as the proposed method.In Experiment 3, we want to show the robustness of the proposed method and the effectiveness of the extracted subtensor features, and then, the classification performance is compared when the subtensor is extracted with different thresholds in the MPCA and MLDA steps.

Experiment 1
Using the same training samples, the classification performance of the proposed method is compared with the NN classification algorithm with vectors and the Wishart algorithm proposed by Lee [1], in which three bands' data are assumed to be statistically independent.The training samples are selected randomly from the ground-truth, and 500 samples are selected for each class.The testing samples are all pixels in the ground-truth.The spatial window size is 3 × 3 in tensor representation.The threshold ρ is set to 0.99 in MPCA, and the threshold is set to 0.995 in MLDA for the proposed algorithm.The classification results are shown in Figure 3. Figure 3a is the classification result of the Wishart algorithm; Figure 3b is the classification result of the Wishart algorithm after the data have been filtered with a 3 × 3 window; Figure 3c is the result of the NN classifier in which the tensor data are unfolded into a vector without dimension reducing; and Figure 3d is the result of the proposed algorithm.The dimension of the original tensor is 9 × 3 × 9.After being processed with MPCA and MLDA, the dimension of the subtensor is 8×3×2.Table 1 lists the classification result, including the correct rates, the total accuracies and the kappa coefficients of the confusion matrices.From Table 1, it can be observed that the classification performance of the proposed algorithm is good.The correct rates of different crops are all higher than 95%, except corn; the overall accuracy is as high as 98.9%; and the kappa coefficient measuring the overall classification effect of the proposed algorithm is 97.9%.The performance of the proposed method is improved for all crops compared with the Wishart classifier and the NN classifier with the original data without dimension reducing, and the accuracies of some crops are improved obviously, such as beet, barley and onions.Although the overall accuracies of the proposed method and that of the Wishart classifier with the filtered data are almost equal, the accuracy rates of the classes of onion and grass are far better than those of the Wishart classifier.It can also be observed that the classification result is better than other algorithms from Figure 3.The burr of the proposed method in homogeneous areas is less than that of the Wishart algorithm.Although the classification accuracies are nearly equal, the classification accuracies are better than the NN classifier without dimension reducing and the Wishart classifier with the filtered data in the areas that are not marked in the ground-truth, especially in the grass area from Figure 3.

Experiment 2
Classification performances of the NN classifier with the two vector learning methods and three tensor learning methods are compared.The two vector learning methods are the PCA and PCA + LDA algorithms; the three tensor learning methods are the MPCA, MLDA and the proposed MPCA + MLDA algorithms.The window size is the same as that used in Experiment 1.To keep the classification performance relatively good, the average dimensions of subtensors of those five algorithms about the same and the values relatively low, the threshold of the PCA step is set to 0.9 in the PCA method, to 0.95 in the PCA+LDA method, the threshold of MPCA only to 0.9, the threshold of MLDA only to 0.99, to 0.97 in the MPCA step and to 0.99 in the MLDA step for the proposed method.The experiment is repeated 100 times.In each trial, the training samples are selected randomly, and the sample number for each crop is 300; the testing samples are all pixels in ground-truth.The average dimensions of the subvectors and subtensors of those algorithms are about the same after being processed, which is 23 for PCA, 32 for MPCA, 21 for MLDA and 21 for the proposed method.The dimension of the PCA + LDA method is determined by the number of classes.Table 2 lists the average classification result of the different algorithms.
According to the results shown in Table 2, firstly, we can observe that the classification performances of the tensor-based learning methods are better than those of the vector-based learning methods.Secondly, for the three tensor-based learning methods, MPCA performs worse than the MLDA and the proposed method when the dimensions of the subtensors are about the same, especially in Types 4 and 11.The performance of the proposed method is slightly better than the MLDA algorithm.Thirdly, the accuracies of the proposed method are improved compared with the Wishart classifier with the filtered data for all crops, except Type 10.

Experiment 3
In the third experiment, the classification performance is compared when the subtensor is extracted with different thresholds in the MPCA and MLDA steps.The threshold is varying from 0.9 to 0.99 for MPCA and varying from 0.92 to 0.999 for MLDA.The training and testing samples are selected the same as in Experiment 2, i.e., 100 trials are carried out for each threshold.The variation of classification accuracies of all 14 crops along with the dimensional change of the subtensors are shown in Figure 4.With the increase of the dimension of subtensors, the accuracies of all crops are improved gradually.The ascent rates are fast when the dimension is less than 20 and tend to be saturated after the dimension exceeds 20, except Crop 10.The variation tendency of the overall accuracy and the kappa coefficient are the same, also increasing rapidly firstly and then growing slowly.It can be observed that the classification accuracy for each crop and the overall accuracy will exceed 90% after the dimension of subtensors is higher than 20, which shows the robustness of the proposed algorithm and the efficiency of the extracted subtensor features.

Conclusions
Taking into consideration the characteristics of frequency, polarization and spatial domain simultaneously for the multi-frequency polarimetric synthetic aperture radar (SAR) data, we have presented a multi-frequency polarimetric SAR classification scheme based on tensor representation and multilinear subspace learning (MSL).Compared with the classic method in which the data or features are represented and learned by vector, the proposed method represents the cell of multi-band data with a third-order tensor in the directions of spatial, frequency and polarization.The principle features are extracted more precisely, and the discrimination ability of the features between classes is improved when using the multi-linear principal component analysis (MPCA) method and the multi-linear extension of linear discriminant analysis (MLDA) method.The proposed algorithm was evaluated with multi-frequency polarimetric SAR data of Airborne Synthetic Aperture Radar (AIRSAR).The classification performances of the three tensor learning algorithms, two vector learning algorithms and the classic Wishart classifier were compared.The accuracies of the proposed algorithm with the dimensional change of subtensors were also tested.Experimental results have demonstrated that the classification accuracy of the proposed algorithm is better than other algorithms; the performance is promising, even when the dimension of the extracted subtensors is low.In addition, the extracted subtensors are three-dimensional, and the overall structure of the original data in the directions of spatial, frequency and polarization is well kept.

= 1 M 1
M m=1 X m .The average tensor of the class c is Xc = Nc Nc i=1 X i , and initialize U I l ×I l .

Figure 2 .
Figure 2. The pseudo-color images and ground-truth: (a-c) Pauli pseudo-color images of C, L and P bands; (d) the ground-truth [21]; (e) the types of crops and the corresponding colors.

Figure 3 .
Figure 3.The classification results: (a) the Wishart algorithm; (b) the Wishart algorithm with filtered data; (c) the NN classifier with original data unfolded into a vector; (d) the proposed algorithm.

Figure 4 .
Figure 4.The accuracy of the proposed algorithm with the dimensional change of the subtensors; (a) correct rates for each crop; (b) overall accuracy and kappa coefficient.

Table 2 .
The classification performance of the PCA only, LDA only, multi-linear principal component analysis (MPCA) only, multi-linear extension of linear discriminant analysis (MLDA) only, the proposed method and Wishart algorithm.