A Hybrid Sparse Representation Model for Image Restoration

Group-based sparse representation (GSR) uses image nonlocal self-similarity (NSS) prior to grouping similar image patches, and then performs sparse representation. However, the traditional GSR model restores the image by training degraded images, which leads to the inevitable over-fitting of the data in the training model, resulting in poor image restoration results. In this paper, we propose a new hybrid sparse representation model (HSR) for image restoration. The proposed HSR model is improved in two aspects. On the one hand, the proposed HSR model exploits the NSS priors of both degraded images and external image datasets, making the model complementary in feature space and the plane. On the other hand, we introduce a joint sparse representation model to make better use of local sparsity and NSS characteristics of the images. This joint model integrates the patch-based sparse representation (PSR) model and GSR model, while retaining the advantages of the GSR model and the PSR model, so that the sparse representation model is unified. Extensive experimental results show that the proposed hybrid model outperforms several existing image recovery algorithms in both objective and subjective evaluations.


Introduction
The purpose of image restoration is to reconstruct high-quality images x from the degraded images y. This is a typical inverse problem, and its mathematical expression is where H denotes the degenerate operator and n is usually assumed to be zero-mean Gaussian white noise. Under different settings, Equation (1) can represent different image processing tasks. When H denotes the identity matrix, Equation (1) represents the image denoising task [1,2]; when H denotes a diagonal matrix with diagonal 1 or 0, Equation (1) represents an image inpainting task [3,4]; when H denotes the blurring operator, Equation (1) represents an image deblurring task [5,6]. In this paper, we focus on the image restoration task.
In order to obtain high-quality reconstructed images, image prior knowledge is usually used to regularize the solution space. In general, image restoration can be expressed as the following minimization problems: where the first term 1 2 y − Hx 2 2 represents data fidelity, the second term R(x) depends on the image prior, and λ is a regularization parameter that balances the two terms. Due to the ill-posed nature of image restoration, the prior knowledge of image plays description, the proposed hybrid sparse representation model is called HSR model. The NSS priors of degraded images are called internal NSS priors, and the NSS priors of external image datasets are called external NSS priors. Figure 1 shows how the HSR model can repair degraded images. The contributions of this paper are summarized as follows: (1) We propose a hybrid sparse representation model that combines the NSS priori of degraded images and external image dataset to make full use of the specific structure of degraded image and the common characteristics of natural image; (2) The introduction of joint model into the HSR not only retains the advantages of the PSR model and GSR model, but also alleviates their respective disadvantages.
The rest of this paper is organized as follows. Section 1 describes the related work of sparse representation. Section 2 introduces how to learn NSS prior from external image corpus. Section 3 introduces the proposed mixed sparse representation model. Section 4 employs an iterative algorithm based on the alternating direction multiplier framework to solve the proposed model. Section 5 presents the experimental results. Section 6 concludes the paper.

Fundamentals of Image Analysis Methods
This section introduces the knowledge of the HSR model. The proposed HSR model uses the NSS prior knowledge of both degraded images and external datasets and introduces a joint model that integrates the PSR model and GSR model. Therefore, the proposed HSR model is based on GSR and PSR. A brief introduction of these two models is given below.

Patch-Based Sparse Representation
The basic unit of patch-based sparse representation (PSR) model is image patch. Given an image x ∈ R N and a dictionary D ∈ R b×M , b ≤ M, where M represents the number of atoms in the dictionary D. The dictionary D in the PSR model is shared.
x i = R i x, ∀i = 1, 2, . . . , n represents the size image patch of √ b × √ b extracted from the position i, and R i represents the extraction operation. The sparse representation of each patch x i is to find the sparse vectors A i with most coefficients zero, that is x i = DA i . The l 0norm represents the number of non-zero elements in a vector. To regularize the parameter matrix A i with l 0 -norm is to expect most of the elements of A i to be 0 and the parameter A i to be sparse. Therefore, by solving the following l 0 -norm minimization problem, each patch x i can be sparsely represented aŝ where 2 denotes l 2 -norm and λ is a regularization parameter. In the image restoration task, the input degraded image y ∈ R b is used because the original image is not available. Extracting image patch y i from degraded images y, each image patch y i can be sparsely represented asÂ In this way, the whole image can be sparsely represented by a set of sparse codes Â i n i=1 .

Group-Based Sparse Representation
Compared with typical PSR models, the GSR model uses patch group as the basic unit of image processing and can produce more promising results in various image processing tasks [3,21,25]. In this subsection, we briefly introduce the GSR model.
Firstly, the image x is divided into n overlapped patch of size . . , n. For each exemplar patch x i , the m most similar matching patches are selected from the search window of size is W × W by the k-nearest neighbor (KNN) method to form the set S G i . Then, all patches in S G i are stacked into a matrix X i ∈ R b×m , with each patch in the collection S G i as a column of the matrix X i , that is is a matrix of all image blocks with similar structures, it is called a patch group, where x i,j represents the j-th similar patch in the i-th patch group. Finally, given a dictionary D i ∈ R b×K , which D is usually learned from each image group, then each patch group X i can be sparsely represented aŝ where B i represents the group sparsity coefficient of each image group, 0 represents the l 0 -norm, and calculate the non-zero items of each column in B i . In image restoration tasks, since the original image is not available, we can only use the input degraded image y ∈ R b . According to the above steps, the image patch y i are extracted from the degraded image y, search for similar matching patches to generate an image group Y i ∈ R b×m , i.e., Y i = {y i,1 , y i,2 , . . . , y i,m }.
The entire image can be sparsely represented by groups of sparse codes B i n i=1 . In the above introduction, y i in the PSR model and Y i in the GSR model are extracted from the same degraded image y.

Learning NSS Priors from External Image Datasets
As mentioned earlier, the traditional sparse representation model only uses the NSS prior of degraded image and ignores the NSS prior of external dataset. In this section, we use the group-based Gaussian mixture model (GMM) [26,30] to learn the external NSS prior from the patch group of a given training image dataset. The following briefly introduces how to learn NSS priors from external image data sets.

Gaussian Component Learning Based on Group GMM
Similar to the construction process of patch groups in Section 2.2, patch groups are extracted from the given external training image dataset, and each patch group is expressed as where e i,j represents the j-th non-local similar patch of the i-th patch group E i , ∀j = 1, 2, . . . , d. In this paper, the GMM model is used to learn k Gaussian components {N(u k , ∑ k )} from the patch group {E i } S i=1 of the external image dataset, and all patches in each patch group are required to belong to the same Gaussian component. The likelihood of a patch group {E i } S i=1 can be expressed as where K is the total number of Gaussian components, µ k is the mean value, ∑ k is the covariance matrix, π k is the weight of Gaussian components, and ∑ K k=1 π k = 1. The GMM model is parameterized by mean vectors {µ k }, covariance matrices {∑ k }, and the weights of Gaussian components {π k }. To facilitate representation, we introduce variables = {µ k , ∑ k , π k } K k=1 . Assume that all patch groups are independent, and the overall objective likelihood function is L = ∏ S i=1 P(E). Taking the log of it, to maximize the objective function of using group-based GMM learning, We can optimize Equation (9) by using the expectation maximization (EM) algorithm [30][31][32]. In the E-step, the posterior probability of the k component calculated by Bayesian formula is In the M-step, for each patch group E i , we update the model parameters as follows By iteratively alternating between the E-Step and M-Step, the model parameters are iteratively updated until convergence is achieved.

Gaussian Component Selection
For patch group Y i of degraded image y, we can select the most appropriate Gaussian component from the training GMM. According to [31], assuming that the image is broken by a Gaussian white noise with a variance of σ 2 e , the covariance matrix of the k-th Gaussian component will be expressed as Σ k + σ 2 e I, where I is the unit matrix. The k-th Gaussian component belonging to the image group Y i can be selected by posterior probability By maximizing Equation (15), the k-th Gaussian component with the highest probability can be selected for each group Y i . Each group E i has the same Gaussian distribution. The covariance matrix of the k-th Gaussian component is denoted by ∑ k . By using the eigenvalue factorization to ∑ k , we have: (16) where U k represents the orthogonal matrix composed of the eigenvector ∑ k and the diagonal matrix Λ k of the eigenvalues. Based on the above GMM learning, the feature vector of U k can represent the statistical structure of NSS changes of natural images, so U k can be used to represent the structural changes of image groups in this component [33]. Finally, we select the best matched U k for each patch group Y i . Since solving the l 0 -norm minimization problem is an NP-hard problem, the l 0 -norm minimization in Equation (6) is replaced by a non-convex l 1 -norm. The sparse model based on external NSS can be expressed as followŝ where C i represents the sparse coefficient of the i-th image group Y i , and ω represents a non-negative constant. After obtaining all the sparse codes Ĉ i n i=1 , a high-quality reconstructed imagex can be obtained.

The Proposed Hybrid Sparse Representation Model
As mentioned above, the traditional sparse representation model only uses the internal NSS priors of degraded images, which leads to over-fitting in the image restoration process. Therefore, this paper uses both the internal NSS priors of degraded images and the external NSS priors of external image dataset. At the same time, the PSR model usually produces some undesirable visual artifacts, and the GSR model leads to over-smoothing effects in various image processing tasks. In order to overcome their shortcomings and improve the image restoration effect, we have introduced a joint model [8] based on both internal and external NSS priors. This model integrates the PSR model and GSR model, instead of using Equations (4) and (6) separately. Combining Equations (4), (6) and (17), the proposed new hybrid sparse representation model is expressed as where N i represents the internal sparse coefficient of the joint sparse representation model, and L i represents the internal joint dictionary of the joint sparse representation model. U k represents the external dictionary, which is learned from the image group of the external image data set using the external NSS prior [26,30]. µ and η represent non-zero constants and act as balance factors to make the solution of Equation (18) more feasible. τ = λ 2 , ϕ = ρ 2 , ω represents the regularization parameter, which is used to balance the sparse coefficients terms of A i 0 , B i 0 , and C i 1 . The sparse coefficient A i corresponds to the sparseness of the image patch on the basis of maintaining the local consistency of the image, which reduces the over-smoothing effect. The sparse coefficient B i corresponds to the sparseness of patch group on the basis of maintaining the non-local consistency of the image and suppresses the undesirable visual artifacts. For specific details of the joint sparse representation model, please refer to [8]. Based on the above analysis, the proposed hybrid sparse representation model not only uses the internal and external NSS priors, but also unifies the sparse representation model.
The hybrid sparse representation model is used in the task of image restoration, and the joint Equations (1) and (18) are expressed as In Equation (19), L represents the internal dictionary of the joint sparse representation model, and U represents the external dictionary. N represents the sparse coefficient of the joint sparse representation model, and C represents the external sparse coefficient. The hybrid sparse representation model proposed in Equation (19) not only comprehensively considers the NSS priors of internal image and external image database, which can provide mutually complementary information for image reconstruction, but also unifies the sparse representation model by combining the PSR model and GSR model.

The Solution Process of the Proposed Hybrid Sparse Representation Model
In this Section, in order to make the proposed model manageable and robust, the alternating direction method of multipliers (ADMM) [34,35] is adopted to solve the largescale optimization problem in Equation (19). Specifically, the minimization of Equation (19) involves three sub-problems, including A i , B i , and C i . Different from the traditional optimization strategies that only considers the fixed values of parameters µ, η, τ, ϕ, and ω, we adaptively adjust all parameters in Equation (19) at each iteration to ensure the stability and practicability of the algorithm. The specific implementation details of the hybrid sparse representation model are given below.

Solution of Hybrid Sparse Representation Model Based on ADMM
Equation (19) is a large-scale non-convex optimization problem. In order to make the optimization problem easy to handle, the alternating direction multiplier method (ADMM) is used. The basic principle of ADMM is to decompose the unconstrained minimization problem into different constrained sub-problems. The following is a brief introduction to ADMM algorithm through a constraint optimization problem, The basic ADMM is shown in Algorithm 1, where t represents the number of iterations.

Algorithm 1. ADMM
Require: Z and N 1: Going back to the hybrid sparse representation model, we transform Equation (19) into two constraint problems and call the ADMM method to solve it. We first transform Equation (19) into an equivalent constraint form by introducing auxiliary variables Z and Q, To facilitate the solution, Equation (21) can be decomposed into two constrained optimization problems, Equation (22) represents the constrained optimization problem of solving the internal joint sparse representation model, and Equation (23) represents the constrained optimization problem of solving the external sparse representation model.

Solution of Internal Sparse Representation Model
Solving the internal sparse representation model in Equation (22) where D represents the fixed dictionary in the PSR model, and D G represents the cascade of all sub-dictionaries D i in the GSR model. Using line 4 in Algorithm 1, we obtain The minimization problem N in Equation (25) is decomposed into A i and B i , and solved respectively, aŝ Finally, using line 5 of Algorithm 1 to update J t , In summary, the minimization of Equation (22) involves three minimization problems, including Z, A i , and B i . Equation (26) represents the PSR model, and Equation (27) represents the GSR model. The implementation details of an effective solution to each sub-problem are given below.

Solution of Sub-Problem Z
Given A and B G , the sub-problem Z in Equation (24) is transformed into, Equation (29) is a quadratic form that has a closed-form solution so that In Equation (30), I denotes the unit matrix of the desired dimension, and J i is the corresponding element in J. In Equation (30), (26) and (27) are used in combination to estimate each.

Solution of Sub-Problem A i
For each image patch in Equation (26), the sub-problem can be re-expressed as where (31) is a sparse representation problem, where the constraint form is directly solved, where θ represents a small constant. Equation (32) can be effectively solved by the orthogonal matching pursuit (OMP) algorithm [36].

Solution of Sub-Problem B i
Given Z and A, the sub-problem in Equation (27) can be transformed into, where R i = Z i − DA i − J i . Solving Equation (33), we find, An important problem in solving sub-problem B i is the choice of dictionary D i . To adapt the local structure of the image, a dictionary based on principal component analysis (PCA) is learned for each group R i . Due to the orthogonality of dictionary D i , Equation (34) can be rewritten asB where R i = D i γ i . We can solve the closed solution for each B i by soft thresholding [31],

Solution of External Sparse Representation Model
Solving the external sparse representation model in Equation (23), defining f (Q) = 1 2η y − HQ 2 2 ,g(C) = ω∑ n i=1 C i 1 , and using line 3 in Algorithm 1, Using line 4 in Algorithm 1, we obtain Finally, using line 5 of Algorithm 1 to update O t , In summary, the minimization of Equation (23) involves two minimization subproblems, including Q and C i . The solution procedure for Q and C i is similar to that in Section 5.2, and the implementation details of an efficient solution for each sub-problem are given below.

Solution of Sub-Problem Q
Given the internal sparse representation model, Equation (23) translates intô Equation (37) is a quadratic form that has a closed-form solution so that In Equation (38), I denotes the unit matrix of the desired dimension and O i is the corresponding element in O.

Solution of Sub-Problem C i
Given Q, the sub-problem C i in Equation (40) can be transformed so that Evidently, Equation (42) can be viewed as a sparse representation problem for each image group Y i . According to Section 2, we can select the best-matching Gaussian component for each group through Equation (15), and then assign the best matching PCA-based dictionary to each group according to Equation (16). Due to the orthogonality of dictionary U i , Equation (42) can be rewritten aŝ where Y i = U i κ i . We can solve the closed solution for each a by soft threshold [37], A complete description of the hybrid sparse representation model for image restoration is given in Algorithm 2.

Algorithm 2. A hybrid sparse representation model for image restoration
Require degraded image y, mask H and group-based GMM 1: Set parameters t, b, W, m, µ, η, τ, ϕ, ω, v,  Select the k-th optimal Gaussian component by Equation (15)  8: Select the dictionary U k by Equation (16)

Adaptive Parameter Adjustment Strategy
There are six parameters in Equation (21), namely µ, η, τ, ϕ, ω, and v. A fixed value is usually chosen for each parameter based on experience. However, this makes it difficult to guarantee the stability and effectiveness of the whole algorithm. To address this problem, an adaptive parameter adjustment scheme is proposed to make the proposed algorithm more stable and practical. An iterative regularization strategy [38] is used to update the estimate of the noise variance σ e . The standard deviation of the noise σ e at the t-th iteration is expressed as where t denotes the number of iterations, δ denotes the scale factor controlling the variance estimation, and the scheme has been widely used for image denoising with Gaussian noise variance estimation [30,38]. Therefore, µ t and η t can be expressed as where γ i and κ i denote the estimated standard deviation ofB i andĈ i [39], respectively. ς denotes a constant with a small value to avoid division by zero. In order to make the proposed algorithm more accurate and practical, according to [37], in the t-th iteration, the ADMM balance factor a is set to v t = 1 where c denotes the scale factor.

Experimental Results
In this section, the experimental results of the proposed HSR model and seven comparison methods are given, including the SALSA [40], BPFA [41], GSR [3], JSR-SR [8], GSRC-NLP [9], IR-CNN [42], and IDBP [43] methods. All experiments were carried out on Intel (R) Core (TM) I7-6700 CPU and 3.40 GHz CPU PC under Matlab 2018B environment. The source code for all competing methods is open source, and we use the default parameter settings. The 13 images used for the experimental tests are shown in Figure 2. In order to evaluate the quality of the restored images, an experimental comparative analysis of the restored images was performed from both objective and subjective aspects. For objective evaluation, the peak signal to noise ratio (PSNR) and structural similarity (SSIM) [44] metrics were used for the experimental comparison of the restored images. The PSNR was calculated as shown in Equations (49) and (50), where X and Y denote the original image and the restored image, respectively, and H × W denotes the size of the image. Equation (20) In Equation (51), SSIM measures similarity in terms of luminance l, contrast c, and image structure s. Where u X and u Y denote the mean of the original image X and the restored image Y of size H × W, respectively; σ X and σ Y denote the variance of the original image X and the restored image Y, respectively; and σ XY denotes the covariance of the original image X and the restored image Y. C 1 , C 2 , and C 3 are constants and introducing a constant can avoid the situation where the denominator is 0. The SSIM indicator is closer to human subjective feelings, and its value range is [0, 1]. The larger the value of SSIM, the more similar the two images are, and the better the effect of image restoration.
For color images, this paper only focuses on the restoration of the luminance channel in YCrCb space. In the group-based GMM learning phase, the training patch group used in the experiment is collected from the Kodak photoCD dataset, which includes 24 natural images.

Objective Evaluation
In the image restoration task, the image restoration results are given for four masks, i.e., 80%, 70%, 60%, and 50% of random pixel loss. The parameters of the HSR model used for image restoration are set as follows: the search window W × W is set to 25 × 25, the size of the image patch is set to 8 × 8, the number of similar patches is set to 60, σ e = √ 2, ς = e −14 , and v = 0.2. We compared the proposed HSR model with seven restoration methods, including SALSA [40], BPFA [41], GSR [3], JPG-SR [8], GSRC-NLP [9], IR-CNN [42] and IDBP [43]. Among these seven methods for comparison, SALSA [40], BPFA [41], GSR [3], JPG-SR [8], and GSRC-NLP [9] are based on traditional image restoration algorithms. The GSR [3], JPG-SR [8], and GSRC-NLP [9] methods are image restoration algorithms based on the traditional GSR model, which belong to the same type of model as our proposed HSR model. SALSA [40] and BPFA [41] are not based on GSR. In order to comprehensively evaluate the performance of the proposed model for image restoration, the proposed HSR model was also compared with algorithms based on deep learning [42,43].
The SALSA model [40] proposes an algorithm belonging to the augmented Lagrangian method family to deal with constraint problems. The method solves optimization problems where the optimal regularization parameters are tuned by manual trial and error, which requires considerable time and effort to achieve the optimal value of the method. The BPFA model [41] utilizes a non-parametric Bayesian dictionary learning method for image sparse representation, and uses image patches as the basic unit of sparse representation, which ignores the similarity between image patches. In terms of the average value, the proposed HSR model is 4.74 dB and 6.19 dB higher than SALAS and BPFA methods respectively. The GSR method [3] is a typical representative of the traditional GSR model, and the JPG-SR method [8] and the GSRC-NLP method [9] are both improved methods based on the GSR model. The GSR method, the JPG-SR method and the GSRC-NLP method only utilize the internal NSS prior. However, the HSR model proposed in this paper combines internal and external NSS priors. In terms of the average value, the HSR model proposed in this paper improves 1.47 dB, 1.43 dB, and 1.06 dB over the GSR, JPG-SR, and GSRC-NLP methods respectively. IRCNN [42] and IDBP method [43] are recovery methods based on deep learning, using powerful prior knowledge of deep neural networks. In terms of the average value, the proposed HSR model improves 3.66 dB and 3.01 dB over the IRCNN and IDBP methods, respectively.
As shown in Tables 1-4, the PSNR of the proposed HSR model on images with a pixel loss rate of 80%, 70%, 60% and 50% is higher than that of SALSA, BPFA, GSR, JRG-SR, GSRC-NLP, IR-CNN and IDBP. It can be seen from the statistical SSIM values in Tables 5-8, that the HSR model is better than other methods in most cases. The experimental results in Tables 1-4 and Tables 5-8 prove that the proposed HSR model is effective and gives good restoration results compared to the comparison method.

Subjective Assessment
The visual comparison between the proposed HSR model in this paper and SALSA [40], BPFA [41], GSR [3], JPG-SR [8], GSRC-NLP [9], IR-CNN [42] and IDBP [43] methods after restoration of the image Mickey with pixel missing rate 80% is given in Figure 3. It can be observed from Figure 3 that the SALSA [40] and BPFA [41] methods cannot recover sharp edges and fine details. The GSR [3] method is better in recovering details, but produces an over-smoothing effect. The JPG-SR [8] method can obtain better visual quality than GSR [3] method. However, the objective evaluation results in Tables 1-5 and Tables 5-8 show that although the JPG-SR [8] method has a higher mean value of PSNR than the GSR [3] method in Tables 1-5, in the actual restoration process, the PSNR and SSIM values of some images after restoration are lower than the restoration results of GSR [3] method. The image restoration effect of the JPG-SR [8] method is unstable, and only some of the image restoration results are better than the GSR [3] method. The GSRC-NLP [9] method can obtain similar visual effects as our proposed HSR model, which is not easy to distinguish from the naked eye. However, according to the experimental results in Tables 1-5 and  Tables 5-8, our proposed HSR model has better objective evaluation results. The visual result of our proposed method is also better in recovering details than the results of IR-CNN [42] and IDBP [43]. The visual results in Figure 3 show that our proposed HSR model retains clear edges and details, especially at higher pixel missing rates, and produces the result with the best visual quality.

Running Time
In this section, we present a comparison of the proposed HSR method with other comparison methods in terms of running time. Taking image Butterfly as an example, the running time of all comparison methods is compared when the pixel loss rate is 50%. As can be seen from Table 9, the processing time of HSR method proposed in this paper is 5000.22 s for the image, which is less than 5027.67 s of the GSRC-NLP method. The proposed HSR method utilizes NSS to construct internal and external image groups and needs to learn the corresponding dictionaries, which requires higher computational workload and therefore consumes more time. To reduce processing time in our future work, learning external NSS priors from the external data set will be done in advance in the Kodak photoCD data set. Through one-time learning from Kodak photoCD data set, the external NSS priors are obtained. The priors learned in advance are applied to speed up the proposed HSR method.

Conclusions
In order to improve the repair performance of the traditional GSR model, we propose a new hybrid sparse representation model. The model uses the NSS prior of degraded image and external image data set, so that the model is complementary in the feature space and the plane. And on this basis, we introduced a joint sparse representation model. The joint model integrates the PSR model and the GSR model, while retaining their advantages, overcoming their shortcomings, and unifying the sparse representation model. Experimental results show that the model is comparable to the test method, and it is better than several state-of-the-art image restoration and map methods in both objective and subjective aspects.