Improvement of Target Detection Based on Signal Dependent Noise Reduction for Hyperspectral Image

In the hyperspectral images (HSI) acquired by the new-generation hyperspectral sensors the signal dependent noise is an important limitation to the detection. Therefore, noise reduction is an important preprocessing step to analyze the information in the hyperspectral image. A signal dependent noise cannot be reduced by conventional linear filtering. Therefore, a new method based on Parallel factor analysis (PARAFAC) decomposition is proposed to estimate the noise of hyperspectral remote sensing image. Then, the estimated noise is used for whitening the colored structural noise. By using this transformation, the data noise from new-generation hyperspectral sensors is diminished, thereby allowing a minimization of negative impacts on hyperspectral detection applications.


Introduction
Hyperspectral sensors collect data in hundreds of narrow contiguous spectral bands, providing powerful means to discriminate different materials by their spectral features.Actually, denoising is of great interest for target detection [1] with the underlying principle of targets being distinguishable.As HSIs are normally produced by a series of sensors, the noise mainly comes from two aspects: signal-independent circuity noise (SI) and signaldependent photon noise (SD) [2].The noise level of HSI may vary dramatically from band to band.The noise variance in each band of HSI is not constant, in particular, there exist some bands at which the atmosphere absorbs so much light that the signal received from the surface is unreliable.Although the SD noise has become as dominant as the SI noise in HSI data collected by new-generation hyperspectral sensors due to the improved sensitivity in the electronic components [2,3,4,5,6,7].For this the denoising methods for those two types of noise are not the same.The SI noise term is generally modelled as additive and spatially stationary in each band, but the variance of the noise varies from band to band.That is to say, the level of the noise is dependent on the average amplitude of each band, but spatially stationary in each band.Based on HSI with high spatial resolution often contain a large number of small homogeneous areas, in [9] an automatic algorithm by dividing an image into several small blocks and calculating local means and local standard deviations of these blocks is developed, then estimating the noise by using a histogram statistical algorithm.To approximate better the noise variance of small blocks, [10] utilized data-masking technology by assuming that image textures are generally smoother than noise.However, the aforementioned methods treat separately each band in HSI data and fail to take spectral information into account.To utilize better the information from high spectral resolution of HSIs, a spectral and spatial decorrelation algorithm [11] based on the multiple linear regression (MLR) model was proposed to estimate noise.The residuals of the MLR model are considered to be noise, while the signal of a pixel at a particular band can be described as a linear combination of the neighboring pixels in the same band and the same spatial pixels in immediately adjacent bands.Based on the mixed noise assumption, the splitting of noise and the original signal from a HSI is usually the first step in these algorithms.Then, maximum likelihood estimation (MLE) is used for the estimation of the SD and SI noise parameters in the second step [12], this method is referred to as the hyperspectral noise parameter estimator (HYNPE) method.In [2], the SD and SI noise were estimated from a single scanning window.Because of the high computational complexity, the algorithm is not widely applied.In summary, MLE-based algorithms suffer from two disadvantages: sensitivity in initial value selection and high computational complexity for the optimal solution.Recently, in [13] denoising algorithm employing a spectral-spatial adaptive total variation model (SSATV), in which the spectral noise differences and spatial information differences are both considered in the process of noise reduction is proposed.In this paper, we propose a novel algorithm that can estimate noise from HSIs with different noise types.
According to the different statistical properties of SI and SD noise, in this paper, we propose a new method based on PARAFAC decomposition [14,15] to remove these two types of noise respectively.Our proposed method is applied to the simulated HSIs in order to evaluate their performances in a controlled environment.Results obtained on the real-world HYperspectral Digital Imagery Collection Experiment (HYDICE) and airborne visible/infrared imaging spectrometer (AVIRIS) Indian Pines HSIs are also presented and discussed.To compare the denoising performance of PARAFAC to other methods, multiway Wiener filter (MWF), prewhiteningmultiway Wiener filter (PMWF) method proposed to reduce colored SI noise [16], HYNPE method, SSATV method and two well-known 2D denoising methods, minimum noise fraction (MNF) and noise-adjusted principal components analysis (PCA) [17], are used in the experiment section.The experimental results show that the proposed method have potential prospective in the reduction of both SI and SD noise in HSIs.
The remainder of the paper is organized as follows : Section 2 reviews some multilinear algebra tools.Section 3 gives the data model of HSI distorted by both SD and SI noise.Section 4 presents the detailed description of our proposed based on multilinear algebra decomposition for enhancement of SNR.Some denoising and comparative results are contained in Section 5. Finally, section 6 concludes the paper.Within the scope of this paper, scalar is denoted by x, vector by x, matrix by X and tensor by §.

Multilinear algebra tools and signal model
In this paper, the order of a tensor is defined as the number of dimensions or modes.Let I 1 , I 2 , I 3 ∈ N denote index upper bound the tensor.A tensor X ∈ R I1×I2×I3 is a real 3-dimensional array, whose element is noted as x i1,i2,i3 , where In the following, some basic multilinear algebra tools used in tensor decompositions are introduced.

Rank-one Tensor
An N -mode tensor X ∈ R I1×I2×I3 being rank 1 means that it can be written as the outer product [18] of 3 vectors, that is: X = a (1) • a (2) • a (3) .So, each element of i2 and a (3) i3 are the i 1 th, i 2 th and i 3 th element of a (1) , a (2) and a (3) , respectively.

n-mode Unfolding
The n-mode vectors are the I n -dimensional vectors obtained from a tensor by varying index i n while keeping the other indices fixed.The so-called n-mode flattened matrix X n ∈ R In×Mn (n = 1, 2, 3) denotes the n-mode unfolding matrix of a tensor X ∈ R I1×I2×I3 , with size I n ×M n where M n = I p ×I q with p = q = n (p, q = 1, 2, 3).The columns of X n are the I n -dimensional vectors obtained from X by varying index i n while keeping the other indices fixed.

n-mode Product
The n-mode product " × n " is defined as the product between a data tensor X ∈ R I1×...×I N and a matrix B ∈ R J×In in mode n.It leads to the tensor U=X × n B of size

Data model
A noisy HSI can be expressed as a third order tensor R ∈ R I1×I2×I3 composed of a multidimensional signal X ∈ R I1×I2×I3 impaired by an additive random noise N ∈ R I1×I2×I3 .The tensor R can be expressed as [19]: where N accounts for both SI and SD noise [20] and its variance depends on the pixel x i1,i2,i3 in the useful signal X.Element-wise, the data model is: where u i1,i2,i3 is a stationary, zero-mean uncorrelated random process independent of x i1,i2,i3 with variance σ 2 u,i3 and w i1,i2,i3 is electronics noise which is zero-mean white Gaussian noise in each band with variance σ 2 w,i3 .The additive term x 1/2 u is the generalized SD noise and denoted as photon noise, w is the SI noise component and is generally assumed to be Gaussian distribution in each band.Then, we can define: and Equation ( 3) can be correspondingly rewritten as With the assumption that x, u and w are independent and both u and w are zero mean and are stationary, the variance of noise N in band i 3 of the HSI could be written as [2,12] where can be expressed as : where X 3 is the 3-mode unfolding matrix of the multidimensional signal tensor X and with N SD,3 and W 3 being the 3-mode unfolding matrices of N SD and W respectively.Using the mean noise variance of the i 3 th spectral band defined as 1/(I 1 I 2 ) is the mean of all x i1,i2,i3 in the i 3 th band of X with i 3 = 1, • • • , I 3 , and the assumption of the independence of x and u, where u is zero-mean and independent between spectral bands, the covariance matrix of the 3-mode unfolding matrix N SD,3 can be expressed as:

Proposed denoising method based on PARAFAC decomposition model
Our aim is to estimate the signal tensor X from tensor R, in the sense of minimum mean square error.From Equation (1) the rank-K s PARAFAC approximation of noisy R is: With the assumption that the noisy tensor R can be exactly expressed by sum of K (K > K s ) rank-1 tensors, then: where k is a residual tensor.Since the rank-1 tensors are orthogonal, then the square error of PARAFAC decomposition According to the definition in Equation (2), a , then a i,k a (2) i,k a (2) so, a Therefore, the minimum of the square error R − R a 2 corresponds to throw away other smaller terms from K s + 1 to K of PARAFAC decomposition.The signal components in the smallest terms from K s + 1 to K are the smallest ones among all the signal components from 1 to K of PARAFAC decomposition.Therefore, for fixed K s , PARAFAC decomposition can reduce the noise, that is to say: the rank-K s PARAFAC approximation of a noisy tensor results in an estimation of the signal,i.e., Since the denoising by PARAFAC decomposition is based on skipping smaller terms from K s + 1 to K where noise components exist, so PARAFAC decomposition has the effect of the reduction of noise, that is to say the residual parts of PARAFAC decomposition, i.e., is an estimate of the noise N .

Estimation of the Optimal Rank K s of PARAFAC Decomposition
In this paper, we assume that the statistical properties of both the signal and the noises are relatively constant in the tensor R and the variance of the SD noise in Equation ( 7) can still be expressed by σ2 u,i3 • μi3 , where σ2 u,i3 and μi3 are the variance of random process u i1,i2,i3 and the mean of all pixels in the i 3 th band of R, respectively.Therefore, the covariance matrix C N can be considered approaching a diagonal matrix.This criterion is used to estimate the optimal rank K s of PARAFAC decomposition for the reduction of noise.

Summary of the Proposed Method for Signal-to-Noise Enhancement
The complete proposed method to reduce both SD and SI noise noise in HSIs can be summarized as follows: 1. Estimate the noise N from R using PARAFAC: (a) Unfold the estimated noise tensor N to n-mode unfolding matrix with n = 1, 2, 3, (b) Calculate the covariance matrix of the n-mode unfolding matrix of the estimated noise tensor, (c) Use the criterion above to assess the result matrix obtained from step b) diagonal matrix.Then, the rank K s of PARAFAC decomposition for the reduction of noise could be obtained from these three steps.2. Estimate the signal tensor X = R − N

Experimental results
To demonstrate the efficacy of different methods for noise reduction, the peak signal to noise ratio (PSNR) of the denoised HSI will be compared.The PSNR can be calculated by PSNR = 10 log 10 (max{X }) 2 MSE (dB) (16) where In addition, the variance of the residual noise, X −X , in the simulated HSI and the variance of the removed noise, R − X , in the real data are compared in this paper after denoising by different considered methods: PARAFAC, PMWF, HYNPE, MWF, SSATV, and two well-known 2D denoising methods, MNF and noise-adjusted PCA.

Simulated Data Experiments
To verify the performance of the proposed method a synthetic HSI is generated according to the data model in Equation ( 3), with the spectral signatures presented in Fig. 1 (a), having size 150 × 150 × 148 .There are six target types and three different spatial sizes 7 × 7 pixels, 2 × 2 pixels and 1 × 1 pixel of each type, which are shown in Fig. 1(b).These targets are mixed to the background by using the linear mixing model with target abundance being 85% (mixing ratio).The random noise is generated with a variance depending on the value of the useful signal according to Equation ( 7) and added into the signal X as Eq. ( 4) to create the simulated HSI data R.In the following sections, the variance of the residual noise in the simulated HSIs after denoising, the peak signal to noise ratio (PSNR) and the ACE target detection results of simulated HSI denoised by the considered methods will be illustrated and discussed.), thus the noise variance of N in this case varies strongly from band to band.The variance of the residual noise, X -X , in the denoised HSI is evaluated at each band and some results are shown in Fig. 2 (b) where the noise variance of N in the raw simulated HSI is also illustrated as a comparison.From Fig. 2 (b) it can be seen that all the considered methods can effectively remove the noise in the simulated HSI since the residual noise variance are much lower than that in raw HSI.But, Fig. 2 (b) demonstrates that the denoised HSI by PARAFAC method contains the least noise when compared against other considered methods.

Target Detection Performance: Probability of Detection
The main purpose of HSI denoising is to improve the results of detection, classification, etc.In this section, we focus on the improvement of target detection using the adaptive coherence/cosine estimator (ACE) [22] which is largely applied to HSI data.The results of ACE target detection of both simulated and real data denoised by  3), the ACE detector can be expressed as: where r j is the vector in the unfolding matrix R 3 of tensor R with j = 1, • • • , I 1 I 2 , s is the target spectrum template Γ is the covariance matrix of R 3 .To assess the performance of detection, the probability of detection (PD) is defined as: and the probability of false alarm (PFA) is defined as: PFA = ns i N fd i ns i (I1×I2−Ni) , where n s is the number of spectral signatures, N i the number of pixels with spectral signature i, N rd i the number of rightly detected pixels, and N fd i the number of falsely detected pixels.To set the values of the two noise variances σ 2 w,i3 , we define some quantities: The signal to noise ratio (SNR) of the synthetic noisy HSI SNR = 10 log 10 X 2 N 2 (dB) (20) then the noise variance of N is σ 2 N = P • 10 − SNR 10 with P = 1/(I 1 I 2 I 3 ) Assuming for the given values of SNR the values of the noise model parameters can be obtained: With the definitions above, it is clear that for δ = 1, σ 2 N SD = σ 2 w , i.e., SD noise source contributes similarly as white SI noise source to the simulated HSI, and there is not a dominant noise source in the simulated HSI, for δ < 1, σ 2 w is higher than σ 2 N SD and the simulated HSI is distorted mainly by white SI noise, otherwise for δ > 1 the SD noise source is dominant.So, in this paper, we consider the cases of δ = [0.1,0.3, 0.5, 1.0, 1.5, 2.0, 2.5] with constant SNR = 30dB which can help to evaluate the influence of SD and SI noise variance on the denoising performance of different methods.Fig. 3 shows the PD of ACE target detection results under the condition PFA=10 −4 for simulated HSI (Fig. 2 (a)) denoised by the considered methods in this paper.In Fig. 3 it is clear that the probability of detection of denoised HSI by PARAFAC method for different values of δ outperforms other methods.The SD and SI noise are removed by throwing away other smallest terms from K s + 1 to K, thus the quality of the denoised HSI by PARAFAC is ameliorated so much that the ACE target detection of the denoised HSI is improved greatly.
For the real-world HSIs, the denosing performance of the proposed method is also verified and discussed in the next section.

Results on Real-World Data
As the previous test is not realistic in the sense that true HSI is simulated according to the data model in Equation ( 3)-( 7), a real-world image is considered in this section.Referred to as HYDICE HSI, was acquired by HYperspectral Digital Imagery Collection Experiment (HYDICE).

Removed Noise Variance and Detection Results
The real-world HYDICE HSI shown in Fig. 4(a) has 150 rows and 140 columns and 148 spectral channels out of 210 with 0.75 m spatial and 10 nm spectral resolution.It can be represented as a 3D data cube, denoted by R ∈ R 150×140×148 .Six targets are added into the HSI and each row of targets in Fig. 4(a) has the same target spectral signature (spectral reflectance) illustrated in Fig. 4(b), which is taken from the image itself.The target size is 5 × 5 pixels along the first column, 3 × 3 pixels along the second one and 1 × 1 pixel along the last one.To assess the denoising results obtained by the proposed method, the removed noise (R − X ) variance is calculated at each band and plotted in Fig. 5 and the receiver operating characteristic (ROC) curves of ACE target detection is presented in Fig. 6.   6 shows the comparison of ROC curves obtained by ACE detector from different methods.It is obvious that the PD values of ACE target detection of denoised HYDICE HSI by PARAFAC method are improved significantly when compared against raw HYDICE HSI counterpart and it is indicated that they are superior to the other considered methods in the reduction of both SD and SI noise.The improvement of target detection of the HYDICE HSI denoised by PARAFAC method is superior to that by SSATV.Since the PARAFAC method has not limitation in denoising colored noise, it can be concluded that the SD noise is at least as dominant as the SI noise in this HSI.The PMWF, MWF, noise-adjusted PCA and MNF methods are not designed for SD noise reduction, so their denoising performance are not ideal, which is reflected indirectly by the target detection results in these ROC curves in Fig. 6.HYNPE method has limitation in removing all SD noise components as shown in Fig. 5, the removed noise by this method is least and correspondingly its denoising performance confines its contribution in the improvement of target detection of this HYDICE HSI.

Conclusion
In this paper, we developed a new multidimensional denoising method based on multiple linear regression and multilinear algebra tools to enhance SNR of HSI data collected by new-generation hyperspectral sensors, distorted by both SI colored and SD noise.To reduce both SD and SI noise, we propose a tensor-based method.PARAFAC decomposition is applied to the noisy HSI to estimate the signal tensor.PARAFAC decomposition must be conducted at the appropriate rank which can be estimated according to the statistical properties of noise.The performance of the proposed PARAFAC method are validated on the simulated HSIs distorted by both SD and SI noise and on the real-world HYDICE.The HYDICE dataset is used to evaluate the denoising and detection results.The experimental results show the efficiency of the proposed denoising algorithm to improve the SNR in hyperspectral images and the target detection.
From the analysis and the comparative study against other similar methods in the experiments, it can be concluded that PARAFAC method can effectively reduce both SD and white or colored SI noise from HSIs.It is also necessary to take into account the noise signal-dependency hypothesis when dealing with HYDICE HSI data.

N
of noise is a diagonal matrix.If the squared norm of the covariance C (3) N 2 is quite close to the sum of the squared diagonal elements I3 i3=1 c 2 i3,i3 , then this C

Figure 1 :
Figure 1: (a) Spectral signatures of the simulated targets and background, (b) Simulated HSI without noise, from top to bottom the index of targets is 1 to 6 respectively

Figure 2 :
Figure 2: (a) Noisy image with SNR =30 dB, (b) Noise variance of raw simulated HSI and residual noise variances of denoised HSI

Figure 3 :
Figure 3: Probability of detection of simulated HSI denoised by different methods with SNR= 30dB and for PFA= 10 −4

Figure 4 :Figure 5 :
Figure 4: (a) HYDICE HSI used to compare different denoising methods, from top to bottom the index of targets is 1 to 6 respectively.(b) Spectral signatures (spectral reflectances) of the simulated targets.(c) Noise variance of raw HYDICE HSI

Figure 6 :
Figure 6: ROC curves of the denoised HYDICE HSIs obtained by ACE detector