Self-Supervised Hyperspectral Inpainting with the Optimisation inspired
  Deep Neural Network Prior
Hyperspectral Image (HSI)s cover hundreds or thousands of narrow spectral
bands, conveying a wealth of spatial and spectral information. However, due to
the instrumental errors and the atmospheric changes, the HSI obtained in
practice are often contaminated by noise and dead pixels(lines), resulting in
missing information that may severely compromise the subsequent applications.
We introduce here a novel HSI missing pixel prediction algorithm, called Low
Rank and Sparsity Constraint Plug-and-Play (LRS-PnP). It is shown that LRS-PnP
is able to predict missing pixels and bands even when all spectral bands of the
image are missing. The proposed LRS-PnP algorithm is further extended to a
self-supervised model by combining the LRS-PnP with the Deep Image Prior (DIP),
called LRS-PnP-DIP. In a series of experiments with real data, It is shown that
the LRS-PnP-DIP either achieves state-of-the-art inpainting performance
compared to other learning-based methods, or outperforms them.
Shuo Li, Mehrdad Yaghoobi
# Self-Supervised Hyperspectral Inpainting with the Optimisation inspired Deep Neural Network Prior

###### Abstract

Hyperspectral Image (HSI)s cover hundreds or thousands of narrow spectral bands, conveying a wealth of spatial and spectral information. However, due to the instrumental errors and the atmospheric changes, the HSI obtained in practice are often contaminated by noise and dead pixels(lines), resulting in missing information that may severely compromise the subsequent applications. We introduce here a novel HSI missing pixel prediction algorithm, called Low Rank and Sparsity Constraint Plug-and-Play (LRS-PnP). It is shown that LRS-PnP is able to predict missing pixels and bands even when all spectral bands of the image are missing. The proposed LRS-PnP algorithm is further extended to a self-supervised model by combining the LRS-PnP with the Deep Image Prior (DIP), called LRS-PnP-DIP. In a series of experiments with real data, It is shown that the LRS-PnP-DIP either achieves state-of-the-art inpainting performance compared to other learning-based methods, or outperforms them.

 Hyperspectral Imaging, Image Inpainting, Plug and Play Denoiser, Self-Supervised Learning.

## I Introduction

In hyperspectral imagery, sensor failures and malfunctions of the HSI acquisition system may result in missing pixels/lines or some spectral bands, significantly hindering the subsequent processing of observed HSI [1]. Hyperspectral inpainting is the task of filling in the missing areas with plausible contents. However, the inpainting of HSIs is a more challenging task than RGB images as each pixel to be filled in is a complex vector with rich spatio-spectral information. Traditional methods such as [2, 3, 4] and [5] either fail when the whole spectral bands of pixels are missing, or their performance is severely compromised if there are a large number of missing pixels. In this work, we treat hyperspectral inpainting as a special case of the reconstruction problem, whose objective is to recover the ground truth from the degraded/masked incomplete images. The low rankness and sparsity of the underlying clean HSI are used here as the priors during reconstruction. The sparse representation (SR) and low rankness (LR) priors have been successfully applied in a wide range of hyperspectral imaging applications such as classification [6], denoising [7], and un-mixing [8]. Recently, researchers have found that the missing spectrum of the HSIs can be predicted through learning on a large dataset [9] or learning from the image itself [10]. The latter is an extension of the Deep Image Prior (DIP) [11] applied to HSIs, achieving the state-of-the-art performance. In [11], authors reveal that the structure of a generative network is sufficient to capture plenty of low-level image statistics prior to any learning. The well-designed CNN networks have been tested on a wide range of tasks such as image denoising, super-resolution and image inpainting, showing promising results which are even competitive to the state-of-the-art deep models trained on large datasets. The "free of external training data" property of the DIP, makes it well-suited for HSI inpainting. The very recent works in [12, 13] raised the point that some trained or untrained neural networks can be directly plugged into the iterative solver to achieve better reconstruction accuracy. Keeping this Plug-and-Play (PnP) idea in mind, one may design a better HSI inpainting algorithm by taking advantage of both the traditional and deep learning techniques. Considering the high computational cost of end-to-end training on extensive HSI data, it is here proposed to use DIP, a self-supervised framework that is free of external training data.

### _Contribution_

This paper aims to develop an effective HSI inpainting algorithm that enjoys the specific learning capability of deep networks, called inductive bias, but does not need any external training data, i.e. self-supervised learning. Specific contributions are the following:

* A novel self-supervised HSI inpainting method, called Low Rank and Sparsity Constraint Plug and play (LRS-PnP), which is presented to solve the most challenging scenarios where the whole spectral bands are missing.
* The use of deep neural networks in replacing the rank constrained-optimisation problem, showing the potential of DIP in learning intrinsic low-rank characteristics of the HSIs.
* A deep hyperspectral prior-based model which better exploits the intrinsic characteristics of HSI data, for achieving state-of-the-art performance.
* Extensive experiments on real data to verify the superiority of proposed LRS-PnP and LRS-PnP-DIP algorithms over existing inpainting solutions.

The rest of paper is organized as follows: section II introduces the proposed framework, section III provides the implementation details and experimental setups, in section IV, the results are discussed. Finally, section V concludes the paper.
## II Proposed Method

The HSI inpainting task can be interpreted as the reconstruction of the clean image \(X\) from its noisy and incomplete measurement \(Y\), in the presence of additive noise \(N\) and masking operator \(M\):

\[Y=M\{X\}+N \tag{1}\]

The clean image \(X\in\mathbb{R}^{q}(q=n_{r}\times n_{c}\times n_{b})\), where \(n_{r},n_{c}\) are the spatial dimension of image, and \(n_{b}\) stands for the total numbers of spectral bands. The operator \(M:\mathbb{R}^{q}\rightarrow\mathbb{R}^{q}\) is the binary mask, with 0 representing the missing pixel and 1 representing the observed and valid pixel. Therefore \(M\) can be presented with a diagonal matrix that has only one per row, which we note it by \(\mathrm{M}\). \(N\) is the additive Gaussian noise of appropriate size. Specifically, \(M\) is often given. The formulation (1) is a linear system which can be written as:

\[\mathbf{y}=\mathrm{M}\mathbf{x}+\mathbf{n} \tag{2}\]

Where \(\mathbf{x},\mathbf{y},\mathbf{n}\) are the vectorized forms of \(X,Y,N\), respectively, and \(\mathrm{M}\) is a diagonal matrix. We first introduce an operator \(P_{i}(\mathbf{x})\) that extracts each i-th patch from the image \(\mathbf{x}\). \(P_{i}(\mathbf{x})\) may cover only the valid pixels or may include the missing pixels depending on the size of \(P_{i}(\cdot)\). We apply the sparse representation on each image patch \(P_{i}(\mathbf{x})\). The inpainted image \(\mathbf{x}^{\star}\) can be obtained by solving the following optimization problem:

\[\begin{split}(\mathbf{x}^{\star},\mathbf{\alpha}^{\star})=& \operatorname*{argmin}_{\mathbf{x},\mathbf{\alpha}}\gamma\|\mathbf{y}-\mathrm{M}\mathbf{x}\|_ {2}^{2}+\mathrm{w}_{\mathrm{lr}}\|\mathbf{x}\|_{*}+\mathrm{w}_{\mathrm{s}}\|\mathbf{ \alpha}\|_{1}\\ &\text{s.t.}\quad\mathbf{x}=\Phi\mathbf{\alpha}\end{split} \tag{3}\]

The first term is the data fidelity term, which we weigh with the parameter \(\gamma\). Due to the ill-posed nature of estimating \(\mathbf{x}\) from \(\mathbf{y}\) only using data fidelity term, the solution is often not unique. For this reason, we introduce another two "priors" to regularize the inpainting problem, namely low rank and sparsity constraints. The second term penalizes the solution \(\mathbf{x}\) to be of low rank, which is often used as the surrogate for the rank minimization problem. The third term constrains the missing pixels to be generated from the subspace approximated by the valid pixels. Similarly, we weigh these two terms with parameters \(w_{lr}\) and \(w_{s}\). The sparse representation problem is solved with a known dictionary \(\Phi\) that is learned only from the noisy and incomplete pixels, or it is a sparsifying transform. By adopting the augmented Lagrangian and introducing the auxiliary variable \(\mathbf{u}\)[14], problem (2) can be rewritten as:

\[\begin{split}(\mathbf{x}^{\star},\mathbf{\alpha}^{\star})=& \operatorname*{argmin}_{\mathbf{x},\mathbf{\alpha}}\gamma\|\mathbf{y}-\mathrm{M}\mathbf{x}\|_ {2}^{2}+\mathrm{w}_{\mathrm{lr}}\|\mathbf{u}\|_{*}+\mathrm{w}_{\mathrm{s}}\sum_{ \mathrm{i}}\|\mathbf{\alpha}_{\mathrm{i}}\|_{1}\\ &+\frac{\mathbf{\mu}_{1}}{2}\|\sum_{i}(P_{i}(\mathbf{x})-\Phi\mathbf{\alpha}_{ i})+\frac{\mathbf{\lambda}_{1}}{\mathbf{\mu}_{1}}\|_{2}^{2}\\ &\text{s.t.}\quad\mathbf{x}=\mathbf{u}\end{split} \tag{4}\]

Where \(\mathbf{\lambda}_{1}\) and \(\mathbf{\mu}_{1}\) are the Lagrangian multiplier and penalty term, respectively. With the help of the alternating direction method of multipliers (ADMM), problem (4) can be solved by the sequential updates of three variables: \(\mathbf{\alpha}\), \(\mathbf{u}\) and \(\mathbf{x}\).

1) _Fixing \(\mathbf{u}\) and \(\mathbf{x}\), and updating \(\mathbf{\alpha}\)_:

\[\begin{split}\mathbf{\alpha}^{k+1}=\operatorname*{argmin}_{\mathbf{\alpha }}\frac{\mathbf{\mu}_{1}^{k}}{2}\sum_{i}\|(P_{i}(\mathbf{x}^{k})+\frac{\mathbf{\lambda}_{1} ^{k}}{\mathbf{\mu}_{1}^{k}})-\Phi\mathbf{\alpha}_{i}\|_{2}^{2}\\ +w_{s}\sum_{i}\|\mathbf{\alpha}_{i}\|_{1}\end{split} \tag{5}\]

which is a patched-based sparse coding problem which can be solved using iterative solvers. In our algorithm, we adopt the PnP-ISTA [15], which has shown promising results over conventional ISTA [16]. Denote the first term in equation (5) as \(f=\frac{\mathbf{\mu}_{1}^{k}}{2}\sum_{i}\|(P_{i}(\mathbf{x}^{k})+\frac{\mathbf{\lambda}_{1 }^{k}}{\mathbf{\mu}_{1}^{k}})-\Phi\mathbf{\alpha}_{i}\|_{2}^{2}\). The whole process can then be replaced by an off-the-shelf denoiser \(\mathcal{D}\) acting on the gradient of \(f\), as it is proposed in [15]. Every single iterate takes the form:

\[\mathbf{\alpha}^{k+1}=\mathcal{D}(I-\nabla f)(\mathbf{\alpha}^{k}) \tag{6}\]

2)_Fixing \(\mathbf{\alpha}\) and \(\mathbf{x}\), and updating \(\mathbf{u}\)_:

\[\mathbf{u}^{k+1}=\operatorname*{argmin}_{\mathbf{u}}w_{lr}\|\mathbf{u}\|_{*}+\frac{\mathbf{ \mu}_{2}^{k}}{2}\|(\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}})- \mathbf{u}\|_{2}^{2} \tag{7}\]

which can be solved by the Singular Value Thresholding (SVT) algorithm. The element-wise soft shrinkage is applied to the singular value of \((\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}})\), as follows,

\[\mathbf{u}^{k+1}=SVT(\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}}) \tag{8}\]

In the proposed LRS-PnP-DIP algorithm, update step of (8) is replaced by a untrained randomized weight neural network \(f_{\theta}(\mathbf{z})\), where \(\theta\) represents the network weights to be updated, and the input \(\mathbf{z}\) is set to be \(\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}}\). i,e, the latent image from the previous iterations:

\[\mathbf{u}^{k+1}=f_{\theta}(\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}}) \tag{9}\]

3) _Fixing \(\alpha\) and \(\mathbf{u}\), and updating \(\mathbf{x}\)_:

\[\begin{split}\mathbf{x}^{k+1}=\operatorname*{argmin}_{\mathbf{x}}\gamma\| \mathbf{y}-\mathrm{M}\mathbf{x}\|_{2}^{2}+\sum_{\mathrm{i}}\|(\mathrm{P}_{\mathrm{i}}( \mathbf{x})+\frac{\mathbf{\lambda}_{1}^{k}}{\mathbf{\mu}_{1}^{k}})-\Phi\mathbf{\alpha}_{1}^{ \mathrm{k}+1}\|_{2}^{2}\\ +\frac{\mathbf{\mu}_{2}^{k}}{2}\|(\mathbf{x}+\frac{\mathbf{\lambda}_{2}^{k}}{ \mathbf{\mu}_{2}^{k}})-\mathbf{u}^{k+1}\|_{2}^{2}\end{split} \tag{10}\]

Closed-form solution for \(\mathbf{x}\) exists as follows:

\[\begin{split}\mathbf{x}^{k+1}=(\gamma\mathrm{M}^{\mathrm{T}}\mathrm{M} +\mathbf{\mu}_{1}^{k}\sum_{\mathrm{i}}\mathrm{P}_{\mathrm{i}}^{\mathrm{T}}\mathrm{P }_{\mathrm{i}}+\mathbf{\mu}_{2}^{k}\mathrm{I})^{-1}\\ (\gamma\mathrm{M}^{\mathrm{T}}\mathbf{y}+\mathbf{\mu}_{1}^{k}\sum_{ \mathrm{i}}\mathrm{P}_{\mathrm{i}}\Phi\mathbf{\alpha}_{1}^{\mathrm{k}+1}+\mathbf{\mu}_{2 }^{k}\mathbf{u}^{k+1}-\sum_{\mathrm{i}}\mathrm{P}_{\mathrm{i}}\mathbf{\lambda}_{1}^{k} -\mathbf{\lambda}_{2}^{k}\mathrm{I})\end{split} \tag{11}\]

4) _Lagrangian and penalty terms updating_:

\[\begin{split}\mathbf{\lambda}_{1}^{k+1}=\mathbf{\lambda}_{1}^{k}+\mathbf{\mu}_{ 1}^{k}(\mathbf{x}^{k+1}-\Phi\mathbf{\alpha}^{k+1})\\ \mathbf{\lambda}_{2}^{k+1}=\mathbf{\lambda}_{2}^{k}+\mathbf{\mu}_{2}^{k}(\mathbf{x} ^{k+1}-\mathbf{u}^{k+1})\end{split} \tag{12}\]
\[\mathbf{\mu}_{1}^{k+1} =\mathbf{\rho}_{1}\mathbf{\mu}_{1}^{k} \tag{13}\] \[\mathbf{\mu}_{2}^{k+1} =\mathbf{\rho}_{2}\mathbf{\mu}_{2}^{k}\]

The proposed Low-Rank and Sparsity Plug-and-Play (LRS-PnP) inpainting model is presented in Algorithm 1.

```
0: masking matrix: \(\mathrm{M}\), noisy and incomplete HSI: \(\mathbf{y}\), learned dictionary: \(\Phi\). denoiser: \(\mathcal{D}\), max iteration: \(It_{max}\).
0:implanted HSI image \(X\).
1:Initialization:\(\mathbf{\lambda}_{1},\mathbf{\lambda}_{2},\mathbf{\mu}_{1},\mathbf{\mu}_{2},\mathbf{\rho}_{1}, \mathbf{\rho}_{2}\).
2:while Not Converged do
3: for \(i=1:It_{max}\) do: \(\mathbf{\alpha}^{k+1}=\mathcal{D}(I-\nabla f)(\mathbf{\alpha}^{k})\)
4:\(\mathbf{u}^{k+1}=SVT(\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}})\)
5: update \(\mathbf{x}\) by ((11)).
6: update Lagrangian parameters and penalty terms.
7:endwhile
```

**Algorithm 1** (LRS-PnP) Algorithm

By replacing the SVT with DIP \(f_{\theta}\), we end up with an extension of the LRS algorithm, denote as LRS-PnP-DIP, is presented in Algorithm 2:

```
0: masking matrix: \(\mathrm{M}\), noisy and incomplete HSI: \(\mathbf{y}\), learned dictionary: \(\Phi\). denoiser: \(\mathcal{D}\), max iteration: \(It_{max}\). DIP: \(f_{\theta}\)
0: inpainted HSI image \(\mathbf{x}\).
1:Initialization DIP parameters, \(\mathbf{\lambda}_{1},\mathbf{\lambda}_{2},\mathbf{\mu}_{1},\mathbf{\mu}_{2},\mathbf{\rho}_{1}, \mathbf{\rho}_{2}\).
2:while Not Converged do
3: for \(i=1:It_{max}\) do: \(\mathbf{\alpha}^{k+1}=\mathcal{D}(I-\nabla f)(\mathbf{\alpha}^{k})\)
4: update \(\theta\) in DIP, with the target \(\mathbf{y}\) and input \(\mathbf{x}^{k}+\frac{\mathbf{\lambda}_{2}^{k}}{\mathbf{\mu}_{2}^{k}}\).
5: update \(\mathbf{x}\) by (11).
6: update Lagrangian parameters and penalty terms.
7:endwhile
```

**Algorithm 2** (LRS-PnP-DIP) Algorithm

## III Implementation Details

We test the proposed inpainting model on the Chikusei airborne hyperspectral dataset, which was taken by Headwall Hyperspec-VNIR-C imaging sensor [17]. The hyperspectral test image has 128 spectral bands, with each 36x36 pixels. The size of the dictionary \(\Phi\) is 1296x2000 which was learned only based on the noisy and incomplete HS image. All images are corrupted with Gaussian noise with the fixed noise strength \(\sigma=0.12\). The mask \(\mathrm{M}\) is applied to all the spectral bands. We use BM3D and Non-local-Mean (NLM) as the plug-and-play denoiser. The implementation of DIP follows the same structures as in Deep Hyperspectral Prior paper [10]. To prevent over-fitting in DIP, we use the early stopping criterion proposed in [18], which automatically detects the near-peak PSNR point using windowed moving variance (WMV). We set \(w_{s}/w_{lr}\) to be 1 and \(\gamma\) to be 0.5. Note that the choice of \(\gamma\) is highly related to the noise level of \(Y\). If the noise level is low, the recovered image \(X\) should be close to noisy observation \(Y\), then parameter \(\gamma\) should be large, and vice versa. For ADMM parameters, \(\lambda_{1}\) and \(\lambda_{2}\) are set to be 0, \(\mu_{1}\), \(\mu_{2}\) are set to be 1, and \(\rho_{1}\), \(\rho_{2}\) are set to be 1, meaning a fixed update step. We use Adam optimizer, and the learning rate is set to be 0.1. We adopt two widely used indicators: Mean Signal-to-Noise Ratio (MPSNR) and Mean Structural Similarity (MSSIM), to evaluate the performance in all experiments.

## IV Experimental Results

The performance of proposed algorithms are compared with the existing traditional methods LRTV [3], and FastHyIn [2], and learning-based method DIP [11], DeepRED [19] and PnP-DIP [20]. The numerical results are reported in Table \(I\) and \(II\). For the DIP, DeepRED, PnP-DIP and the proposed LRS-PnP-DIP, the same U-net in [10] was used as the backbone. For the fairness of comparison, all associated parameters with U-net are kept fixed, including the noise standard deviations and regularization strength. The MPSNR and MSSIM are obtained through running algorithms for 20 times on the same test image. For comparison with FastHyIn and LRTV, 25% of the pixels were randomly masked using a uniform distribution for the missing bands. The inpainting performance is shown in Figure 1.

It shows that LRS-PnP algorithm can generate a more consistent and realistic spectrum in the missing region compared to the learning-based methods. LRS-PnP captures local structures, such as a sudden change in the materials, e.g. see the first image. In contrast, methods such as DeepRED and PnP-DIP tend to generate much smoother contents in the non-missing areas. However, there are severe distortions and artefacts of the missing regions. We believe this is due to the nature of DIP which mainly focuses on learning global features and characteristics. For this reason, it is often used with other regularizers such as Total Variation (TV) to preserve more details [21]. It is worth mentioning that the inpainting result of the double deep prior algorithm LRS-PnP-DIP is visually and qualitatively better than those of other inpainters.

## V Conclusion

The novel hyperspectral inpainting algorithms called LRS-PnP and its extension self-supervised LRS-PnP-DIP which can effectively handle missing pixels from noisy and incomplete HS images in the most challenging scenario where the whole spectrum bands are missing. The new methods exploit spectral and spatial redundancy of HSIs and require no training data except, the test image. A comparison of LRS-PnP and LRS-PnP-DIP with the state-of-the-art algorithms is conducted on a real HS image, leading to the conclusion that LRS-PnP yields similar while LRS-PnP-DIP yields better performance against other learning-based methods, while not using a pre-training
step using a large set of HS images. As the future work, one direction is to accelerate the proposed algorithms to achieve a cost-efficient and real-time HSI inpainting. Another direction of interest is to explore and optimize the training of DIP, especially when it is used in the loop, as mentioned in a lines of DIP-related works [12, 13, 20, 22, 23]. The other direction is to theoretically show the convergence of the algorithm, which we left for a future work.
