# Research on Hyperspectral Image Reconstruction Based on GISMT Compressed Sensing and Interspectral Prediction

Authors: Sheng Cang; Achuan Wang
Journal: International Journal of Optics
 (2020)

---

## Abstract

Hyperspectral remote-sensing images have the characteristics of large transmission data and high propagation requirements, so they are faced with transmission and preservation problems in the process of transmission. In view of this situation, this paper proposes a spectral image reconstruction algorithm based on GISMT compressed sensing and interspectral prediction. Firstly, according to the high spectral correlation of hyperspectral remote-sensing images, the hyperspectral images are grouped according to the band, and a standard band is determined in each group. The standard band in each group is weighted by the GISMT compressed sensing method. Then, a prediction model of the general band in each group is established to realize the remote-sensing image reconstruction in the general band. Finally, the difference between the actual measured value and the predicted value is calculated. According to the prediction algorithm, the corresponding difference vector is obtained and the predicted measured value is iteratively updated by the difference vector until the hyperspectral reconstructed image of the relevant general band is finally reconstructed. It is shown by experiments that this method can effectively improve the reconstruction effect of hyperspectral images.


---

## Body

## 1. Introduction

Hyperspectral images are a set of images with spectral resolutions in the order of10−2λ magnitude, which typically contain hundreds of spectral bands and each of which is imaged separately. In the ultraviolet to mid-infrared region, the hyperspectral image is imaged by the hyperspectral sensor in the satellite for each target region in each band. This kind of acquisition of hyperspectral images can obtain rich spectral information and spatial information in the target area. However, due to the high acquisition band of the hyperspectral spectrum, the data of the hyperspectral image are huge, and the transmission and the processing are difficult. Therefore, how to effectively solve the problems of acquisition, transmission, and preservation of hyperspectral images has become the key to the effective application of hyperspectral technology.Among the existing hyperspectral reconstruction methods, Zhao [1] proposed a method of compressed sensing reconstruction of hyperspectral images based on variable sampling rates, providing a new idea of using different sampling rates for the reference band and the nonreference band [1]. Wand and Fend [2] proposed a multihypothesis prediction compressive sensing reconstruction algorithm based on the space spectrum association. For the reconstruction of nonreference bands, image blocks of nonreference bands were predicted according to the adjacent images of nonreference bands and the reconstructed images of reference bands [2]. In the literature [3], a fast sparse decomposition algorithm based on image decomposition on the overcomplete atomic library was proposed. The greedy clustering algorithm was used to remove the redundancy between the spectra of images and strengthen the correlation of each spectrum band [3]. Jia et al. [4] proposed a prediction method of reference bands to ordinary bands, by adding the prediction of reference bands in adjacent groups to ordinary bands [4]. These methods also play a guiding and helpful role in the study of this paper. However, these methods fail to make full use of the spectral correlation information between various bands of hyperspectral images to improve the quality of reconstruction, and for the reconstruction of reference bands, they cannot effectively reduce the sampling rate of information.Therefore, in view of this situation, this paper designs a hyperspectral image reconstruction method based on GISMT compression sensing and interspectral prediction. The significance of this paper mainly has two aspects: the innovation of the sampling mode and the innovation of image reconstruction [5]. First of all, the images of every band of the hyperspectrum are not reconstructed. In fact, hyperspectral data are grouped according to the spectral correlation of the hyperspectrum. One standard band is determined for each group, while the rest are general bands. Secondly, the GISMT algorithm is used to reconstruct the image data of standard bands to obtain the reconstructed images of each group of standard bands. The GISMT algorithm is a kind of algorithm based on compression perception, compressed sensing theory as a new kind of effective mechanism of image acquisition and sampling method, which is different from traditional. It has the following characteristics: (1) the signal sparse sex is no longer subject to the bandwidth size limit and depends on the sparse signal itself, namely, the content and structure of the important information in the signal; (2) signal sampling and compression processing are completed simultaneously; (3) while retaining complete key information, fewer measurements are obtained after sparse representation; and (4) it can reconstruct the original signal from a small number of measured values with high quality and express the signal structure completely. Compared with the common image compression processing mechanism, the difference lies in that the pure image compression mechanism needs to meet certain conditions and preconditions, and compress and discard redundant data to reduce information quantity and ensure image quality [6]. However, the compression sensing theory does not need to satisfy the criteria of the Nyquist sampling theorem and directly conduct the sampling and compression process of signals. The simple and efficient GISMT algorithm optimizes and solves the improved joint sparse representation model to realize the remote-sensing image reconstruction algorithm [7]. By adjusting the optimal p value, the redundancy between data is reduced to a certain extent, and the number of observations needed for signal reconstruction is greatly reduced. The image reconstruction of the hyperspectral reference band achieved by the GISMT algorithm is simpler and more efficient in calculation, which can not only reduce the computational complexity but also better suppress the noise and ringing effect and protect the image edge details. Then, a spectral prediction model is established, and the prediction of the general bands in each group is carried out according to the model, so as to get the difference between the predicted values and the collected values of the actual general bands [8]. Finally, the difference vector of the difference value is reconstructed, the predicted value is updated according to the vector, and the original image of the band is restored. This algorithm makes full use of the characteristics of hyperspectral correlation, reduces the amount of data acquisition, and through experimental verification, this algorithm can improve the effect of hyperspectral image reconstruction.


## 2. Background

### 2.1. Spatial Correlation of Hyperspectral Images

The hyperspectral image correlation reflects the degree of correlation between adjacent pixels in the spatial position of the image [9]. In this paper, a hyperspectral remote-sensing image is selected and partially magnified 16 times, as shown in Figure 1(a). Subjectively, the gray value between the locally adjacent pixels does not change much. Objectively, the histogram of the gray distribution of the part is shown in Figure 1(b). The abscissa is the gray value of the image, and the ordinate is the frequency at which the grayscale appears. It can be clearly seen from the histogram that the gray value between the local adjacent pixel points is concentrated between 43 and 90, which indicates the hyperspectral remote-sensing image has a certain correlation in the spatial direction.Figure 1
(a) 16 times hyperspectral remote-sensing image. (b) Grayscale values of hyperspectral remote-sensing images.
(a)
(b)


### 2.2. Interspectral Correlation of Hyperspectral Images

The interspectral correlation of hyperspectral data refers to the correlation of adjacent band pixels in the same spatial position, which includes statistical correlation and structural correlation: the pixel values of each band image of hyperspectral data are the same regional features and the reflection value of each band [10]. Because the wavelengths of adjacent bands are similar, the adjacent band images have similar gray values, which is called statistical correlation, and its size depends largely on the spectral resolution. Also, the ground objects have the same goal, and they have the same spatial topology, which is called structural correlation [11]. It mainly analyzes statistical correlation here. Assuming that the pixel values of the band m and the band n in the p-th row and the q-th column are f(p, q, m) and f(p, q, n), respectively, the interspectral correlation coefficient between the band m and the band n is as follows:(1)Rm,n=∑p=1K∑q=1Ffp,q,m−umfp,q,n−un∑p=1K∑q=1Ffp,q,m−um2∑p=1K∑q=1Ffp,q,n−un2,where um and un are the mean values of the gradations of the m-th band and the n-th band image, respectively. In this paper, a hyperspectral remote-sensing image is selected, as shown in Figure 2(a), which is the 45th band grayscale image, and Figure 2(b) is the 49th band grayscale image. Calculated by using formula (1), the correlation coefficient of these two bands of grayscale is 0.99653, so the spectral correlation of the remote-sensing images in the visible light is high, which provides the basic guarantee for the subsequent interspectral prediction.Figure 2
(a) 45th band grayscale image. (b) 49th band grayscale image.
(a)
(b)


### 2.3. Compressed Sensing Theory

The trend of “three features” in remote-sensing data development brings a high data rate to a remote-sensing image [12]. In order to reduce complexity of calculation and increase the speed of calculation, this paper uses the thought of block compressed sensing and TV method to realize the reconstructive algorithm of the remote-sensing image [13]. Firstly, it is based on block sampling, dividing the original remote-sensing image into some nonoverlapped image pieces of the same size, and acquiring an improved joint sparse representation model by basic prior knowledge. Random sampling is implemented in each subimage. Finally, the TV method and ALM reconstructive algorithm are used to get few of measurements to rebuild the original remote-sensing image from all pieces of the image in compressed sampling.An application premise of the compressed sensing theory is that signals must be sparse or compressible. A discrete signal of lengthN is set as follows: xm, m∈1,2,3,…,M according to the signal theory, the signal xm can be represented by a linear combination of basis ψ=ψ1,ψ2,ψ3,…,ψM:(2)x=∑m=1Mψmam=ψa,where x,a∈RM∗1, ψ∈RM∗M. When signal xm has only ψ nonzero coefficients K≪M on the basis ψ, am is called the sparse basis of signal xm and the signal is considered to have sparsity or compressibility. The commonly used sparse bases are sine basis, cosine basis, and so on.In the coding measurement model of compressed sensing, the observed values obtained are not samples of the signal itself, but the projected values of the signal from high dimension to low dimension. The random projection process of the signal is a very critical step in the theory of compressed sensing, which greatly reduces the sampling rate of the signal. Each observation is a combinatorial function of the original signal.Suppose that signalx is projected onto a set of observation vector ϕ=ψ1,ψ2,ψ3,…,ψMT to get the observation value ym=x∣ψmT; that is,(3)y=ϕx,where x∈RM∗1, y∈RN∗1, and ϕ∈RM∗N. Substituting equation (2) into equation (3), we can get(4)y=ϕx=ϕψa=Θa.Among them,Θ=ϕψ=RN∗M. Since the observed dimension is much smaller than the signal dimension, x needs to be solved first by solving a in equation (4) and then substituting it into equation (2) to solve x.The compressed sensing theory points out that if the original signal is to be reconstructed accurately, the measurement matrix must satisfy the principles of incoherence and limiting isovolumency.


## 2.1. Spatial Correlation of Hyperspectral Images

The hyperspectral image correlation reflects the degree of correlation between adjacent pixels in the spatial position of the image [9]. In this paper, a hyperspectral remote-sensing image is selected and partially magnified 16 times, as shown in Figure 1(a). Subjectively, the gray value between the locally adjacent pixels does not change much. Objectively, the histogram of the gray distribution of the part is shown in Figure 1(b). The abscissa is the gray value of the image, and the ordinate is the frequency at which the grayscale appears. It can be clearly seen from the histogram that the gray value between the local adjacent pixel points is concentrated between 43 and 90, which indicates the hyperspectral remote-sensing image has a certain correlation in the spatial direction.Figure 1
(a) 16 times hyperspectral remote-sensing image. (b) Grayscale values of hyperspectral remote-sensing images.
(a)
(b)


## 2.2. Interspectral Correlation of Hyperspectral Images

The interspectral correlation of hyperspectral data refers to the correlation of adjacent band pixels in the same spatial position, which includes statistical correlation and structural correlation: the pixel values of each band image of hyperspectral data are the same regional features and the reflection value of each band [10]. Because the wavelengths of adjacent bands are similar, the adjacent band images have similar gray values, which is called statistical correlation, and its size depends largely on the spectral resolution. Also, the ground objects have the same goal, and they have the same spatial topology, which is called structural correlation [11]. It mainly analyzes statistical correlation here. Assuming that the pixel values of the band m and the band n in the p-th row and the q-th column are f(p, q, m) and f(p, q, n), respectively, the interspectral correlation coefficient between the band m and the band n is as follows:(1)Rm,n=∑p=1K∑q=1Ffp,q,m−umfp,q,n−un∑p=1K∑q=1Ffp,q,m−um2∑p=1K∑q=1Ffp,q,n−un2,where um and un are the mean values of the gradations of the m-th band and the n-th band image, respectively. In this paper, a hyperspectral remote-sensing image is selected, as shown in Figure 2(a), which is the 45th band grayscale image, and Figure 2(b) is the 49th band grayscale image. Calculated by using formula (1), the correlation coefficient of these two bands of grayscale is 0.99653, so the spectral correlation of the remote-sensing images in the visible light is high, which provides the basic guarantee for the subsequent interspectral prediction.Figure 2
(a) 45th band grayscale image. (b) 49th band grayscale image.
(a)
(b)


## 2.3. Compressed Sensing Theory

The trend of “three features” in remote-sensing data development brings a high data rate to a remote-sensing image [12]. In order to reduce complexity of calculation and increase the speed of calculation, this paper uses the thought of block compressed sensing and TV method to realize the reconstructive algorithm of the remote-sensing image [13]. Firstly, it is based on block sampling, dividing the original remote-sensing image into some nonoverlapped image pieces of the same size, and acquiring an improved joint sparse representation model by basic prior knowledge. Random sampling is implemented in each subimage. Finally, the TV method and ALM reconstructive algorithm are used to get few of measurements to rebuild the original remote-sensing image from all pieces of the image in compressed sampling.An application premise of the compressed sensing theory is that signals must be sparse or compressible. A discrete signal of lengthN is set as follows: xm, m∈1,2,3,…,M according to the signal theory, the signal xm can be represented by a linear combination of basis ψ=ψ1,ψ2,ψ3,…,ψM:(2)x=∑m=1Mψmam=ψa,where x,a∈RM∗1, ψ∈RM∗M. When signal xm has only ψ nonzero coefficients K≪M on the basis ψ, am is called the sparse basis of signal xm and the signal is considered to have sparsity or compressibility. The commonly used sparse bases are sine basis, cosine basis, and so on.In the coding measurement model of compressed sensing, the observed values obtained are not samples of the signal itself, but the projected values of the signal from high dimension to low dimension. The random projection process of the signal is a very critical step in the theory of compressed sensing, which greatly reduces the sampling rate of the signal. Each observation is a combinatorial function of the original signal.Suppose that signalx is projected onto a set of observation vector ϕ=ψ1,ψ2,ψ3,…,ψMT to get the observation value ym=x∣ψmT; that is,(3)y=ϕx,where x∈RM∗1, y∈RN∗1, and ϕ∈RM∗N. Substituting equation (2) into equation (3), we can get(4)y=ϕx=ϕψa=Θa.Among them,Θ=ϕψ=RN∗M. Since the observed dimension is much smaller than the signal dimension, x needs to be solved first by solving a in equation (4) and then substituting it into equation (2) to solve x.The compressed sensing theory points out that if the original signal is to be reconstructed accurately, the measurement matrix must satisfy the principles of incoherence and limiting isovolumency.


## 3. The Algorithm of Image Reconstruction

### 3.1. Correlation Analysis and Band Selection of Hyperspectral Remote-Sensing Images

The spatial properties of hyperspectral images are very similar to those of ordinary natural images. The spectral images of each band can be regarded as a common grayscale image, but there is a strong correlation between these grayscale images. In this paper, based on the strong interspectral correlation of the hyperspectral image and combined with the reconstruction method of the gray image, a hyperspectral image compression sensing reconstruction method based on the correlation between spectra is proposed:(1)
Calculating the interspectral correlation of hyperspectral images. SupposeY∈QM×M×n be only the hyperspectral experimental image. ym,n,r represents the pixel value of the pixel position (m, n) of the r-th band, and yr¯ represents the average value of the pixel values of all the pixels in the r-th band. The calculation formula is as follows:(5)yr¯=1M×M∑m=1M∑n=1Mym,n,r.(i)
Supposey˜m,n,r=ym,n,r−y˜r, and the calculation formula of the correlation between the r-th band and the r+1-band is as shown in equation (6).(ii)
It is shown thatc∈0,1.(6)cr,r+1=∑m=1M∑n=1My˜m,n,ry˜m,n,r+1∑m=1M∑n=1My˜m,n,r2∑m=1M∑n=1My˜m,n,r+12.(2)
Grouping the bands of hyperspectral images according to the interspectral correlation. First, the total number of groupst is determined according to the total number of bands of the hyperspectral image, and the size of the group is h. There is one standard band and h−1 general bands in each group; that is, a total of t standard bands are determined. The research shows that taking 10 bands is most suitable. If the number of bands is less than 10, the reconstruction is performed according to the actual number of bands. Because of the band correlation of hyperspectral remote-sensing images, the second band of each group is generally selected as the standard band.


### 3.2. Sparse Representation of the Reference Band and Random Projection

The reference band is processed first. The random projection matrix selects the Gaussian measurement matrix commonly used in the theory of compressed sensing. The results show that the Gauss matrix satisfies the constrained isovolumicity and are irrelevant to almost any sparse or compressible signal, fully complying with the requirements of random projection matrices. Suppose the Gaussian measurement matrixϕ∈QM×N, and it obeys normal distribution with a mean of 0 and a variance of 1/M. p1,1,p2,1,…,pt,1 are t reference bands, which are all sampled at a large sampling rate.In this paper, a biorthogonal wavelet basisψ is selected to represent hyperspectral images sparsely. Biorthogonal wavelet basis is a generalization of the concept of orthogonal wavelet basis. In order to overcome the shortcoming of the orthogonal wavelet, the biorthogonal wavelet is introduced to construct a finite-length filter with linear phase. That is, equations (2) and (4) can be reduced as follows:(7)ym=ϕψam,1,m=1,2,…,t,where ymm=1,2,3,…,t is the measured value. For simplified expression, suppose Y=ymm=1,2,…,t, and then equation (4) is reduced as follows:(8)Y=ϕψam,1,m=1,2,…,t.In the random projection of hyperspectral images, the measurement matrix and the measurement method used in each band of each group are the same, so the image domain has a strong spectral correlation hyperspectral image, which is also very high in its measurement domain.


### 3.3. Reference Band Image Reconstruction Algorithm

In the reconstruction process of the standard band, it is only necessary to reconstruct the standard band of each group since each group has a standard band. This paper designs a GISMT reconstruction algorithm based on compressed sensing.Assuming anN-dimensional signal x ∈ RN, the transform coefficient vector of x in the Ψ domain can be expressed as ΨTx. If there are only K nonzero elements (N >> K) or only K larger components in the transform coefficient ΨTx, the signal x can be considered to be K sparse or approximately K sparse. Suppose there is a randomly projected M × N-dimensional (M << N) measurement matrix A. The linear projection of the signal x yields the M-dimensional measurement y, which is mathematically expressed as y = Ax, and the measured value y is the original signal x. The sampled and compressed sampled signals are processed simultaneously, including all the information of the original signal x. The purpose of solving the problem of compressed sensing reconstruction is to reconstruct a complete original signal x by using a small number of measured values y after compression sampling. The mathematical expression is as follows (constrained optimization problem):(9)minxψxps.t.y=Ax,where p is a value of 1 or 0. ∗1 represents a norm which stands for the sum of the absolute values of all elements in the vector. ∗0 represents the 0 norm, which stands for the number of nonzero elements in the vector.Image prior knowledge plays a key role in the performance of image reconstruction algorithms. The remote-sensing image itself has two basic prior knowledge features, namely, local smoothness and nonlocal self-similarity [14].Local smoothness is a key technique for smoothing images directly in the spatial domain, which is based on the premise that the images processed by high-pass filtering are sparse. On this basis, local smoothing processing of image can be achieved. That is, selecting a feature module in the image and finding the similar module according to the similarity of adjacent pixels. The basic models of local smoothing include mean filtering model, median filtering model, Gaussian filtering model, bilateral filtering model, and total variational filtering model. Local smoothness can be represented by 2D.Nonlocal self-similarity was first used in image desiccation. The algorithm uses the similarity between pixels in the neighborhood to remove noise. In the process of processing, the pixel points in a neighborhood can be extended to the corresponding pixel points on different image blocks. Since the image block contains much more information than a single pixel, the nonlocal self-similarity processing of image block is better than that of pixel. Nonlocal self-similarity can also be applied in image reconstruction. That is, for images with high redundancy of image blocks, the image blocks with high redundancy can be superimposed for three-dimensional transformation, that is, 3D representation [15]. The relationship between 2D and 3D is shown in Figure 3.Figure 3
The relationship between 2D and 3D.The local smoothness describes the smoothness of the local region segmentation of the remote-sensing image in the 2D spatial domain, that is, the gray area of the adjacent pixel. There is similarity in the degree values, while the nonlocal self-similarity describes the repeated texture or structural features of the nonlocal regions in the 3D transform domain of the remote-sensing image. In order to realize the high quality reconstruction of remote-sensing images, two effective regularization parameters are designed based on regularization to represent these two basic prior knowledge features. Considering the overall sparsity of remote-sensing images, the expression is more sparse [16]. By controlling the values of the two regularization parameters, respectively, it can more adaptively represent the local 2D sparsity of the spatial domain of the remote-sensing image and the nonlocal 3D sparseness of the transform domain. The improved joint sparse representation model of remote-sensing images is as follows:(10)ψcosx=τ⋅ψL2Dxp+λ⋅ψN3Dxp,where τ and λ indicate the bound term of regularization, which balance the sparsity of two prior knowledge features. The joint sparse is denoted by cosine.Specifically,ψL2Dxp corresponds to local flatness prior knowledge. It maintains the local continuity of the image and can suppress the noise effectively. Suppose the Laplace distribution is presented by the gradient image of the remote-sensing image, and then depict the sparsity of the local 2D spatial domain by utilizing the filtered image with certain sparsity through filtered convolution for the remote-sensing image; the mathematical representation is as follows:(11)ψL2Dxp=Dxp=Dvxp+Dhxp,where D is the gradient operator and Dh and Dv represent horizontal and vertical gradient operators, respectively.Depict the overall self-similarity of nonlocal 3D transform domain through the combination of coefficient sparsity of all image blocks in 3D groups; its the mathematical representation is as follows:(12)ψN3Dxp=Θxp=∑i=1nT3DZxip.In this formula, Θ represents the column vector of all conversion coefficients of imagex. Substitute joint sparse representation formula (10)–(12) into formula (9), and depict two inherent sparsities of the remote-sensing image on the basis of reducing the image sampling rate while retaining characteristics of the original image. Then, the target function based on the reconstruction problem of compressed sensing joint sparse remote-sensing image will be obtained as follows:(13)minxτ⋅Dxp+λ⋅Θxps.t.y=Ax.Based on the regularization thought, the auxiliary variableDx = d will be introduced. Transform the bound optimization problem of formula (10) into unbound optimize problem of model (14) to solve(14)minx,d,u12Ax−y22+α2x−u22+η2Dx−d22+τ⋅dp+λ⋅Θup.In this formula, regularization items areα and η ⟶ ∞, formula (13) has the same solution with formula (14), and alternating minimization will be adopted to find the solution.(1)
x subproblem:D and u are given, and the solution of x can be found using the following formula:(15)minx12Ax−y22+α2x−u22+η2Dx−d22.The approximate solution ofx obtained from expression (13) is as follows:(16)x=F−1FηDTd+FA∗∘FyηFDhTDh+FDvTDv+FA∗∘FA.In the formula,F is two-dimensional Fourier transform, F−1 is inverse two-dimensional Fourier transform, ∗ is conjugate complex numbers, and ∘ is component-specific multiplication.(2)
Subproblem of u:D and x are given, and u can be solved using the following formula:(17)minuα2x−u22+λΘup.Regardingr as observations of some type noises of x, formula (18) can be equivalently transformed as(18)minu12u−r22+λα∑k=1nT3DZukp.In this formula, suppose every element ofu and r follows independent distribution with a kind of zero mean value and σ2 variance; besides, Gaussian process is not needed in the middle. Due to u, r ∈ RN, Θu, and Θr∈RK, the greatest probability of a three-dimensional conversion coefficient vector at each iteration is existed in two equations:(19)1Nuk−rk22=σ2,1KΘuk−Θrk22=σ2.Combining (16) into (15), formula (20) can be acquired:(20)u=minu12Θu−Θr22+KλNαΘup.A general threshold function (GST) is proposed to realize the iterative contraction mechanism, so as to solve the sparse minimization problem of the nonconvex Lp-norm.The general threshold function GST is as follows, and the specific implementation process is given in reference [7]:(21)TpGSTy;λ=0;ify≤τpGSTλ,sgny⋅SpGSTy;λ;ify>τpGSTλ.According to the GST function, the approximate solution of each iterationui is as follows:(22)ui=TpGSTΘri;KλNα.(3)
d subproblem:u and x are given, and define dref = Dx; d can be solved using the following formula:(23)mindη2d−dref22+τdp.Use the GST function, the approximate solution of each iterationdi is as follows:(24)di=TpGSTdiref;τη.In conclusion, a flow chart that summarizes the algorithm is shown in Figure4.Figure 4
Algorithm flowchart based on GISMT remote-sensing image reconstruction.


### 3.4. Nonreference Band Image Reconstruction Algorithm

Image reconstruction of the nonreference band is completed on the basis of image reconstruction of the reference band. First of all, we build a spectral multidirectional prediction model and use this model for prediction of the nonreference band image. The predicted value of the nonreference band can be obtained with the spectral multidirectional prediction model. As a result, we can obtain the difference between the predicted value and the actual measured value, reconstruct the measurement difference, and update and modify the predicted value. Then, iterative operation looms above the reference to the original value of the data, so as to realize for the band of reference image reconstruction.Interspectral multidirectional prediction is used to predict nonreference wave bands in a wave band group. The data source that will be predicted is the reference wave bandpm,1 in this group, the reference wave band pm−1,1 in the last group, and the reference wave pm+1,1 in the next group. This paper predicts the nonreference wave band of the z group through three reconstructed reference wave bands, namely, pm,1, pm−1,1, and pm+1,1; therefore, the prediction model is established as follows:(25)pm,n′=αx˜m,1+βx˜m−1,1+γx˜m+1,1,where α, β, and γ are predictive coefficients, in α=h−f+1/h, β=f−1/2h, and λ=f−1/2h, h is the band number in the group, and f is the band position in the group. Based on the group strategy stated above, the higher the nonreference band is in the group, the more relevant it is to the reference band within the group, in the prediction model, for the nonreference band with the top position, the reference band of the same group has a larger weight; for the nonreference band with lower position, the reference band of the last group and next group has a larger weight. When the first group of bands is used for prediction, β=0; when the last group of bands is used for prediction, γ=0.The specific reconstruction algorithm is as follows.First of all, conduct measure difference, according to the predicted valuepm,n′ solved by formula (21) of the prediction model and utilize the measurement matrix ϕ that is similar to the reference image; the measured value qm,n′=ϕpm,n′ of the predicted value can be obtained, the practical measurement value of the nonreference wave band pi,j obtained through random projection is qm,n=ϕpm,n, and then the difference between two measured values is follows:(26)qm,n′=ϕpm,n−pm,n′=ϕΔp.In the formula,Δp is the difference between the original value pi,j of the current band and the multidirectional predicted value pi,j′, and the corresponding measurement difference is qm,n′, which equals the difference obtained through randomly projecting the difference vector.Afterwards, reconstruct and restore the difference vector from the measured differenceΔq with the GISMT algorithm. The hyperspectral image has a strong interspectral relevance; therefore, after multidirectional interspectral prediction, the predicted value pm,n′ is very close to the original value pm,n; thus, the obtained measured value is more sparse than the original value pm,n, so the reconstructed difference Δq and the restored difference vector Δq are more effective than reconstructing qm,n and restoring pm,n.At last, the measured valuepm,n′ is corrected by reconstructed and restored p^, and then the updated measured value pm,n′=pm,n+Δx^ will be acquired. Obviously, the updated measured value is closer to the original value. Iterate the above stated steps, including the measured difference, reconstruct difference, and correct predicted value, which can allow the predicted value to keep approaching the original value, when the reconstructed difference vector Δp is less than δ; that is, when the predicted value stops updating, cease iteration, and the predicted value at this time can be regarded as a practical value. The reconstructed algorithm of the ordinary image is shown in Figure 5.Figure 5
The reconstructed algorithm of the ordinary image.


## 3.1. Correlation Analysis and Band Selection of Hyperspectral Remote-Sensing Images

The spatial properties of hyperspectral images are very similar to those of ordinary natural images. The spectral images of each band can be regarded as a common grayscale image, but there is a strong correlation between these grayscale images. In this paper, based on the strong interspectral correlation of the hyperspectral image and combined with the reconstruction method of the gray image, a hyperspectral image compression sensing reconstruction method based on the correlation between spectra is proposed:(1)
Calculating the interspectral correlation of hyperspectral images. SupposeY∈QM×M×n be only the hyperspectral experimental image. ym,n,r represents the pixel value of the pixel position (m, n) of the r-th band, and yr¯ represents the average value of the pixel values of all the pixels in the r-th band. The calculation formula is as follows:(5)yr¯=1M×M∑m=1M∑n=1Mym,n,r.(i)
Supposey˜m,n,r=ym,n,r−y˜r, and the calculation formula of the correlation between the r-th band and the r+1-band is as shown in equation (6).(ii)
It is shown thatc∈0,1.(6)cr,r+1=∑m=1M∑n=1My˜m,n,ry˜m,n,r+1∑m=1M∑n=1My˜m,n,r2∑m=1M∑n=1My˜m,n,r+12.(2)
Grouping the bands of hyperspectral images according to the interspectral correlation. First, the total number of groupst is determined according to the total number of bands of the hyperspectral image, and the size of the group is h. There is one standard band and h−1 general bands in each group; that is, a total of t standard bands are determined. The research shows that taking 10 bands is most suitable. If the number of bands is less than 10, the reconstruction is performed according to the actual number of bands. Because of the band correlation of hyperspectral remote-sensing images, the second band of each group is generally selected as the standard band.


## 3.2. Sparse Representation of the Reference Band and Random Projection

The reference band is processed first. The random projection matrix selects the Gaussian measurement matrix commonly used in the theory of compressed sensing. The results show that the Gauss matrix satisfies the constrained isovolumicity and are irrelevant to almost any sparse or compressible signal, fully complying with the requirements of random projection matrices. Suppose the Gaussian measurement matrixϕ∈QM×N, and it obeys normal distribution with a mean of 0 and a variance of 1/M. p1,1,p2,1,…,pt,1 are t reference bands, which are all sampled at a large sampling rate.In this paper, a biorthogonal wavelet basisψ is selected to represent hyperspectral images sparsely. Biorthogonal wavelet basis is a generalization of the concept of orthogonal wavelet basis. In order to overcome the shortcoming of the orthogonal wavelet, the biorthogonal wavelet is introduced to construct a finite-length filter with linear phase. That is, equations (2) and (4) can be reduced as follows:(7)ym=ϕψam,1,m=1,2,…,t,where ymm=1,2,3,…,t is the measured value. For simplified expression, suppose Y=ymm=1,2,…,t, and then equation (4) is reduced as follows:(8)Y=ϕψam,1,m=1,2,…,t.In the random projection of hyperspectral images, the measurement matrix and the measurement method used in each band of each group are the same, so the image domain has a strong spectral correlation hyperspectral image, which is also very high in its measurement domain.


## 3.3. Reference Band Image Reconstruction Algorithm

In the reconstruction process of the standard band, it is only necessary to reconstruct the standard band of each group since each group has a standard band. This paper designs a GISMT reconstruction algorithm based on compressed sensing.Assuming anN-dimensional signal x ∈ RN, the transform coefficient vector of x in the Ψ domain can be expressed as ΨTx. If there are only K nonzero elements (N >> K) or only K larger components in the transform coefficient ΨTx, the signal x can be considered to be K sparse or approximately K sparse. Suppose there is a randomly projected M × N-dimensional (M << N) measurement matrix A. The linear projection of the signal x yields the M-dimensional measurement y, which is mathematically expressed as y = Ax, and the measured value y is the original signal x. The sampled and compressed sampled signals are processed simultaneously, including all the information of the original signal x. The purpose of solving the problem of compressed sensing reconstruction is to reconstruct a complete original signal x by using a small number of measured values y after compression sampling. The mathematical expression is as follows (constrained optimization problem):(9)minxψxps.t.y=Ax,where p is a value of 1 or 0. ∗1 represents a norm which stands for the sum of the absolute values of all elements in the vector. ∗0 represents the 0 norm, which stands for the number of nonzero elements in the vector.Image prior knowledge plays a key role in the performance of image reconstruction algorithms. The remote-sensing image itself has two basic prior knowledge features, namely, local smoothness and nonlocal self-similarity [14].Local smoothness is a key technique for smoothing images directly in the spatial domain, which is based on the premise that the images processed by high-pass filtering are sparse. On this basis, local smoothing processing of image can be achieved. That is, selecting a feature module in the image and finding the similar module according to the similarity of adjacent pixels. The basic models of local smoothing include mean filtering model, median filtering model, Gaussian filtering model, bilateral filtering model, and total variational filtering model. Local smoothness can be represented by 2D.Nonlocal self-similarity was first used in image desiccation. The algorithm uses the similarity between pixels in the neighborhood to remove noise. In the process of processing, the pixel points in a neighborhood can be extended to the corresponding pixel points on different image blocks. Since the image block contains much more information than a single pixel, the nonlocal self-similarity processing of image block is better than that of pixel. Nonlocal self-similarity can also be applied in image reconstruction. That is, for images with high redundancy of image blocks, the image blocks with high redundancy can be superimposed for three-dimensional transformation, that is, 3D representation [15]. The relationship between 2D and 3D is shown in Figure 3.Figure 3
The relationship between 2D and 3D.The local smoothness describes the smoothness of the local region segmentation of the remote-sensing image in the 2D spatial domain, that is, the gray area of the adjacent pixel. There is similarity in the degree values, while the nonlocal self-similarity describes the repeated texture or structural features of the nonlocal regions in the 3D transform domain of the remote-sensing image. In order to realize the high quality reconstruction of remote-sensing images, two effective regularization parameters are designed based on regularization to represent these two basic prior knowledge features. Considering the overall sparsity of remote-sensing images, the expression is more sparse [16]. By controlling the values of the two regularization parameters, respectively, it can more adaptively represent the local 2D sparsity of the spatial domain of the remote-sensing image and the nonlocal 3D sparseness of the transform domain. The improved joint sparse representation model of remote-sensing images is as follows:(10)ψcosx=τ⋅ψL2Dxp+λ⋅ψN3Dxp,where τ and λ indicate the bound term of regularization, which balance the sparsity of two prior knowledge features. The joint sparse is denoted by cosine.Specifically,ψL2Dxp corresponds to local flatness prior knowledge. It maintains the local continuity of the image and can suppress the noise effectively. Suppose the Laplace distribution is presented by the gradient image of the remote-sensing image, and then depict the sparsity of the local 2D spatial domain by utilizing the filtered image with certain sparsity through filtered convolution for the remote-sensing image; the mathematical representation is as follows:(11)ψL2Dxp=Dxp=Dvxp+Dhxp,where D is the gradient operator and Dh and Dv represent horizontal and vertical gradient operators, respectively.Depict the overall self-similarity of nonlocal 3D transform domain through the combination of coefficient sparsity of all image blocks in 3D groups; its the mathematical representation is as follows:(12)ψN3Dxp=Θxp=∑i=1nT3DZxip.In this formula, Θ represents the column vector of all conversion coefficients of imagex. Substitute joint sparse representation formula (10)–(12) into formula (9), and depict two inherent sparsities of the remote-sensing image on the basis of reducing the image sampling rate while retaining characteristics of the original image. Then, the target function based on the reconstruction problem of compressed sensing joint sparse remote-sensing image will be obtained as follows:(13)minxτ⋅Dxp+λ⋅Θxps.t.y=Ax.Based on the regularization thought, the auxiliary variableDx = d will be introduced. Transform the bound optimization problem of formula (10) into unbound optimize problem of model (14) to solve(14)minx,d,u12Ax−y22+α2x−u22+η2Dx−d22+τ⋅dp+λ⋅Θup.In this formula, regularization items areα and η ⟶ ∞, formula (13) has the same solution with formula (14), and alternating minimization will be adopted to find the solution.(1)
x subproblem:D and u are given, and the solution of x can be found using the following formula:(15)minx12Ax−y22+α2x−u22+η2Dx−d22.The approximate solution ofx obtained from expression (13) is as follows:(16)x=F−1FηDTd+FA∗∘FyηFDhTDh+FDvTDv+FA∗∘FA.In the formula,F is two-dimensional Fourier transform, F−1 is inverse two-dimensional Fourier transform, ∗ is conjugate complex numbers, and ∘ is component-specific multiplication.(2)
Subproblem of u:D and x are given, and u can be solved using the following formula:(17)minuα2x−u22+λΘup.Regardingr as observations of some type noises of x, formula (18) can be equivalently transformed as(18)minu12u−r22+λα∑k=1nT3DZukp.In this formula, suppose every element ofu and r follows independent distribution with a kind of zero mean value and σ2 variance; besides, Gaussian process is not needed in the middle. Due to u, r ∈ RN, Θu, and Θr∈RK, the greatest probability of a three-dimensional conversion coefficient vector at each iteration is existed in two equations:(19)1Nuk−rk22=σ2,1KΘuk−Θrk22=σ2.Combining (16) into (15), formula (20) can be acquired:(20)u=minu12Θu−Θr22+KλNαΘup.A general threshold function (GST) is proposed to realize the iterative contraction mechanism, so as to solve the sparse minimization problem of the nonconvex Lp-norm.The general threshold function GST is as follows, and the specific implementation process is given in reference [7]:(21)TpGSTy;λ=0;ify≤τpGSTλ,sgny⋅SpGSTy;λ;ify>τpGSTλ.According to the GST function, the approximate solution of each iterationui is as follows:(22)ui=TpGSTΘri;KλNα.(3)
d subproblem:u and x are given, and define dref = Dx; d can be solved using the following formula:(23)mindη2d−dref22+τdp.Use the GST function, the approximate solution of each iterationdi is as follows:(24)di=TpGSTdiref;τη.In conclusion, a flow chart that summarizes the algorithm is shown in Figure4.Figure 4
Algorithm flowchart based on GISMT remote-sensing image reconstruction.


## 3.4. Nonreference Band Image Reconstruction Algorithm

Image reconstruction of the nonreference band is completed on the basis of image reconstruction of the reference band. First of all, we build a spectral multidirectional prediction model and use this model for prediction of the nonreference band image. The predicted value of the nonreference band can be obtained with the spectral multidirectional prediction model. As a result, we can obtain the difference between the predicted value and the actual measured value, reconstruct the measurement difference, and update and modify the predicted value. Then, iterative operation looms above the reference to the original value of the data, so as to realize for the band of reference image reconstruction.Interspectral multidirectional prediction is used to predict nonreference wave bands in a wave band group. The data source that will be predicted is the reference wave bandpm,1 in this group, the reference wave band pm−1,1 in the last group, and the reference wave pm+1,1 in the next group. This paper predicts the nonreference wave band of the z group through three reconstructed reference wave bands, namely, pm,1, pm−1,1, and pm+1,1; therefore, the prediction model is established as follows:(25)pm,n′=αx˜m,1+βx˜m−1,1+γx˜m+1,1,where α, β, and γ are predictive coefficients, in α=h−f+1/h, β=f−1/2h, and λ=f−1/2h, h is the band number in the group, and f is the band position in the group. Based on the group strategy stated above, the higher the nonreference band is in the group, the more relevant it is to the reference band within the group, in the prediction model, for the nonreference band with the top position, the reference band of the same group has a larger weight; for the nonreference band with lower position, the reference band of the last group and next group has a larger weight. When the first group of bands is used for prediction, β=0; when the last group of bands is used for prediction, γ=0.The specific reconstruction algorithm is as follows.First of all, conduct measure difference, according to the predicted valuepm,n′ solved by formula (21) of the prediction model and utilize the measurement matrix ϕ that is similar to the reference image; the measured value qm,n′=ϕpm,n′ of the predicted value can be obtained, the practical measurement value of the nonreference wave band pi,j obtained through random projection is qm,n=ϕpm,n, and then the difference between two measured values is follows:(26)qm,n′=ϕpm,n−pm,n′=ϕΔp.In the formula,Δp is the difference between the original value pi,j of the current band and the multidirectional predicted value pi,j′, and the corresponding measurement difference is qm,n′, which equals the difference obtained through randomly projecting the difference vector.Afterwards, reconstruct and restore the difference vector from the measured differenceΔq with the GISMT algorithm. The hyperspectral image has a strong interspectral relevance; therefore, after multidirectional interspectral prediction, the predicted value pm,n′ is very close to the original value pm,n; thus, the obtained measured value is more sparse than the original value pm,n, so the reconstructed difference Δq and the restored difference vector Δq are more effective than reconstructing qm,n and restoring pm,n.At last, the measured valuepm,n′ is corrected by reconstructed and restored p^, and then the updated measured value pm,n′=pm,n+Δx^ will be acquired. Obviously, the updated measured value is closer to the original value. Iterate the above stated steps, including the measured difference, reconstruct difference, and correct predicted value, which can allow the predicted value to keep approaching the original value, when the reconstructed difference vector Δp is less than δ; that is, when the predicted value stops updating, cease iteration, and the predicted value at this time can be regarded as a practical value. The reconstructed algorithm of the ordinary image is shown in Figure 5.Figure 5
The reconstructed algorithm of the ordinary image.


## 4. Analysis of Experiment and Its Results

The experimental data are HIS data of HJ-1A star launched in September 2008, this paper chooses the hyperspectral remote-sensing data of HIS 2 class products in the Da Xing’an Ling of Heilongjiang province as the data source, and the band number was 255. This paper chooses the hyperspectral image numbered 3070541 to perform the experiment. Demonstrating the second image and the third image of the seventh group of No. 1019906 image, namely, the 72nd and 73rd bands, among these bands, the 72nd band is the standard band of the seventh group, while the 73rd band is the ordinary band of this group. The measurement matrix of this experiment is the Gaussian measurement matrix, the sampling rate of the reference band is 0.21, while the ordinary band is 0.07, and the sparse matrix is a biorthogonal wavelet base. This paper chooses the MT-Bayes algorithm that has a good reconstructed effect at present and the GPSR algorithm as a comparison experiment, and the grouped reference band and the nonreference band method will also be used for reconstruction (the specific implementation process of GPSR is presented in [4]). This paper will present the reconstructed effect of the 72nd and the 73rd band of the original image, and the 72nd reconstructed effect is shown in Figure 5. The reconstructed effect of the 73rd band is shown in Figure 6. From Figure 6, it can be perceived that due to the high sampling rate of the reference image, the reconstructed effect that the three methods exerts on the reference image has a slight difference. It can be seen from Figure 7 that the reconstructed effect of interspectral multidirectional prediction used in this paper is much better than the former two algorithms in reconstructing nonreference images.Figure 6
The comparison of the reconstructed effect of the 72rd reference band. (a) Original. (b) MT-Bayes. (c) GPSR. (d) Proposed algorithm.
(a)
(b)
(c)
(d)Figure 7
The comparison of the reconstructed effect of the 73rd reference band. (a) Original. (b) MT-Bayes. (c) GPSR. (d) Proposed algorithm.
(a)
(b)
(c)
(d)As shown in Figure6, the original image has high resolution and prominent details. Figure 6(b) is close to the original image, and the reconstruction effect is good. Figure 6(c) has a vague outline, while Figure 6(d) has few details. Visually, this algorithm is superior to the other two algorithms.As shown in Figure7, the original image has high resolution and clear contour. Figure 7(b) is close to the original drawing with perfect details. Figure 7(c) is not clear in detail, and Figure 7(d) is partially distorted. Visually, the proposed algorithm has more details and clear outlines.The reconstructed effect can be measured by PSNR (peak signal-to-noise ratio), and the calculation formula is as follows:(27)PSNR=20log101MSE.In this formula, MSE is mean square error, and the calculation formula of MSE is as follows:(28)MSE=1mn∑i=1m∑j=1nxi,j−x^i,j2.In this formula,m and n are the width and height of the single band of hyperspectral image, respectively, xi,j is the pixel value that the original image at point (i, j), and x^i,j is the pixel value that the reconstructed image at point (i, j).This paper chooses a band group to calculate PSNR with the 67th band as the reference image. The comparison of PSNR results of the reconstructed image between the MT-Bayes algorithm, GPSR algorithm, and proposed algorithm is shown in Table1.Table 1
The effect comparison of PSNR.
Band image
Algorithm
72
73
74
75
76
77
78
MT-Bayes
33.58
27.36
26.69
26.44
27.81
26.51
27.34
GPSR
33.71
28.19
28.58
29.08
28.86
27.71
28.23
Proposed algorithm
34.34
32.89
33.47
33.21
34.17
33.20
33.76It can be perceived from Table1 that the algorithm put forward in this paper has a better reconstructed effect in the 72nd band of the standard band. From the PSNR of other nonreference bands in this group, it can be seen that, with low sampling rate, the interspectral multidirectional prediction model presented in this paper prominently promotes the reconstructed effect of the nonreference image.The selection of the sampling rate of the reference and nonreference bands has a great influence on the performance of the reconstructed algorithm. We define the average sampling rateR=1/N×rh+N−1/N×rl, where N is the size of the band group, rh is the sampling rate of the reference band, and rl is the sampling rate of the nonreference band. This paper compared the average PSNR of three reconstructed algorithms by selecting different average sampling rates, and the results are shown in Table 2.Table 2
The effect comparison of PSNR.
0.12
0.17
0.22
0.27
0.32
0.37
MT-Bayes
28.84
30.36
32.16
33.88
35.77
37.71
GPSR
29.35
31.26
32.95
34.35
36.21
37.91
Proposed algorithm
33.88
35.70
37.11
35.33
39.38
42.17It can be found from Table2 that compared with the GPSR algorithm and the MT-Bayes algorithm, the average PSNR improved 4-5 dB under the different average sampling rates with the grouping and prediction methods for the reference image presented in this paper; therefore, the reconstructed effect of the algorithm proposed in this paper under different sampling rates raises prominently.


## 5. Conclusion

Certain difficulties have been brought to the transmission process due to large data volume and more band data of the hyperspectral image. To this end, this paper proposed a new hyperspectral reconstructed algorithm based on GISMT compressed sensing and interspectral prediction. First, the group collected hyperspectral remote-sensing image data according to the interspectral relevance features of the hyperspectral image and confirmed a standard band in each group, and other bands are ordinary bands. The group adopted a new GISMT algorithm based on compressed sensing to reconstruct the image type. Then, they established an interspectral prediction model in the light of interspectral relevance. This model mainly predicts the difference between the target bands and the former and latter bands. After that, they reconstructed a remote-sensing image of the ordinary band in each group according to this model. The experiment proved that compared with reconstructed algorithms of other hyperspectral remote-sensing image, this algorithm can remarkably improve the reconstructed effect of the image, which has certain practical application significance.



---

*Source: 7160390-2020-01-29.xml*