袁宇麗，呂俊瑞，羅學(xué)剛,3,4

（1.內(nèi)江師范學(xué)院 計算機科學(xué)學(xué)院，四川 內(nèi)江 641100；2.攀枝花學(xué)院 數(shù)學(xué)與計算機學(xué)院，四川 攀枝花 617000；3.四川旅游學(xué)院 信息與工程學(xué)院，四川 成都 610100；4.電子科技大學(xué) 信息與通信工程學(xué)院，四川 成都 610054）

引言

條帶噪聲是星載和機載多探條帶測器光譜儀圖像中常見的結(jié)構(gòu)噪聲。高光譜成像傳感器在成像過程中不可避免地會對采集到的HSI 數(shù)據(jù)產(chǎn)生條帶噪聲。HSI 中的條帶噪聲將嚴(yán)重影響高光譜圖像的視覺質(zhì)量，對混解、圖像融合、目標(biāo)檢測、分類等后續(xù)處理提出了很大的挑戰(zhàn)[1]。為了獲得更好的視覺質(zhì)量，提高HSI 后續(xù)應(yīng)用的能力，去除HSI 中的條帶噪聲是一個重要研究課題。針對去條帶的問題，學(xué)者們提出了眾多模型[1-5]。其中，基于全變分模型[6-8]的圖像分離方法實現(xiàn)HSI 條帶噪聲去除效果良好，受到學(xué)者們的廣泛關(guān)注，但在實際應(yīng)用中依然存在結(jié)構(gòu)細(xì)節(jié)丟失等不足。為了增強去噪適應(yīng)性，大量改進的方法先后被提出。2016 年Bouali 等人提出了一種自適應(yīng)單向變分的HSI 條帶去除（unidirectional total variation, UTV）[7]模型，該方法利用方向TV 刻畫條帶結(jié)構(gòu)，分離條帶效果良好，但對于無條帶區(qū)域的識別不夠準(zhǔn)確，特別是對嚴(yán)重條帶噪聲的識別性能較差；Dou 等人在2018 年提出基于方向?0范數(shù)稀疏約束的去條帶噪聲模型[9]，采用方向?0范數(shù)的正則化來描述條帶全局稀疏分布，但全局稀疏度不能反映條帶分量的內(nèi)部結(jié)構(gòu)，當(dāng)條帶稠密時導(dǎo)致條帶稀疏性失效；2019 年Wang等人提出RBSD（ reweighted block sparsity destriping, RBSD）去條帶模型[10]，利用加權(quán)塊稀疏和UTV 正則化約束提高圖像變化區(qū)域的局部平滑度，由于該模型沒有考慮高光譜影像的光譜特征，性能仍有較大的提升空間；Liu 等人[11]在2019 年提出基于小波域的低秩組稀疏去條帶模型（ wavelet-domain low-rank/group-sparse destriping,WLRGSD），在變換域下的低秩稀疏約束分離條帶，效果良好，但條帶和圖像結(jié)構(gòu)在小波域中都處于高頻區(qū)域，不能較好地辨識容易導(dǎo)致圖像模糊；2020 年Yang 等人[12]提出聯(lián)合單向變分和非凸低秩構(gòu)建（TV and nonconvex low-rank, TVNLR）去條帶模型，利用高階變分較好地自適應(yīng)圖像結(jié)構(gòu)，但沒有考慮條帶分量結(jié)構(gòu)特點和高光譜影像光譜特征，容易引起非周期條帶分離不完全的缺陷。

為了避免條帶噪聲去除過程中丟失影像細(xì)節(jié)，本文提出一種基于加權(quán)塊稀疏聯(lián)合最小最大非凸懲罰約束的HSI 條帶噪聲去除方法（WBS-MCP）。WBS-MCP 在保持影像邊緣和加強區(qū)域平滑性方面具有良好的優(yōu)勢，特別是對于非周期條帶噪聲去除獲得了更好的效果。

1 WBS-MCP 模型及優(yōu)化求解

本節(jié)首先描述條帶噪聲模型，然后介紹MCP函數(shù)和WBS-MCP 模型及優(yōu)化求解過程，最后描述算法流程。

1.1 條帶噪聲描述

設(shè)定遙感高光譜圖像為Y∈RM×N×B，空間大小為M×N（M為寬度，N為高度），光譜波段數(shù)為B，將HSI 條帶效應(yīng)建模為加性噪聲，模型建立如（1）式：

其中：Y表示含噪HSI 圖像；X表示無噪聲的HSI 圖像；S表示條帶噪聲；Y、X和S具有相同維數(shù)，大小為M×N×B。HSI 條帶噪聲去除就是從觀測噪聲圖像中還原出具有圖像本身特征潛在影像X。為了更好地描述條帶噪聲的低秩和稀疏性，以圖像塊為單位將HSI 數(shù)據(jù)從沿條帶垂直方向和光譜方向利用滑動窗口分解構(gòu)建圖像矩陣，描述形式如（2）式：

式中Y、X和S分別表示HSI 影像按圖像塊分組排列的觀測圖像矩陣、干凈圖像矩陣和條帶噪聲矩陣。

1.2 MCP 懲罰

UTV[7]對圖像梯度利用?1范數(shù)進行約束，體現(xiàn)條帶噪聲的稀疏性特征。研究發(fā)現(xiàn)?0范數(shù)為最佳描述稀疏性，然而由于?0范數(shù)求解困難，采用?1范數(shù)和核范數(shù)等凸函數(shù)近似描述低秩和稀疏，獲得近似?0范數(shù)的效果。在圖像恢復(fù)中，非凸懲罰通常比凸近似逼近具有更好的性能，能更加準(zhǔn)確描述低秩和稀疏。因此，非凸懲罰函數(shù)廣泛應(yīng)用在描述低秩和稀疏的應(yīng)用中，其中常見的有bridge[13]、capped-?1[14]、平滑剪切絕對偏差（SCAD）[15]和極小極大非凸懲罰（MCP）[16]。文獻[17]對上述約束函數(shù)在圖像恢復(fù)中進行數(shù)值實驗比較，發(fā)現(xiàn)MCP 模型在效率和精度上均優(yōu)于經(jīng)典的TV 模型，性能比其他幾種懲罰函數(shù)表現(xiàn)得更突出[17]。

對于一維變量t, MCP 標(biāo)量函數(shù)定義如（3）式[17]：

式中 α＞0 和 β≥1分別為兩個常數(shù)。

擴展到二維矩陣K∈Rp×q， MCP 函數(shù)定義為其中σi(K)表示矩陣K的SVD 奇異值,r表示秩。

那么，相應(yīng)的MCP 閾值操作算子可表示為

圖1 為?0、?1和MCP 3 種懲罰函數(shù)的閾值示意圖。從圖中可以發(fā)現(xiàn)，MCP 函數(shù)具有連續(xù)性、稀疏性和無偏性等特征，保證MCP 正則化問題具有良好的統(tǒng)計性質(zhì), 表示非凸稀疏恢復(fù)MCP 模型比凸稀疏恢復(fù)模型具有更好的性能。由于MCP 函數(shù)具有良好的性質(zhì)，MCP 函數(shù)已經(jīng)在穩(wěn)健主成分分析問題上成功應(yīng)用。

圖1 3 種懲罰函數(shù)的閾值示意圖Fig.1 Threshold schematic of 3 penalty functions

1.3 WBS-MCP 模型描述

本節(jié)描述一種將MCP 約束和加權(quán)?2,1范數(shù)融入到UTV 模型中，聯(lián)合空間-光譜下條帶低秩性和稀疏性正則化項的去條帶新模型WBS-MCP：

其中：λ1、λ2、λ3為非負(fù)參數(shù)；||S||w,2,1表示加權(quán)?2,1范數(shù)，定義為和分別表示水平和垂直方向的微分算子。（5）式的第1 項和第2 項分別描述條帶噪聲低秩性和稀疏性約束，第3 項為估計圖像水平方向的TV 正則化項。圖2 詳細(xì)地描述了WBS-MCP 去條帶模型。

1.4 模型優(yōu)化求解

WBS-MCP 模型將（5）式分解為多個子問題采用ADMM 算法迭代求解，通過引入兩個輔助變量Q=ΔyS,D=ΔxX，（5）式可以重寫為

圖2 WBS-MCP 去條帶模型Fig.2 Stripe removal model of WBS-MCP

由于模型中目標(biāo)函數(shù)的4 個變量具有可分性，采用增廣拉格朗日函數(shù)可重寫為

其中：G1、G2分別為與約束的拉格朗日乘數(shù)； ρ為懲罰參數(shù)。

將WBS-MCP 模型的復(fù)雜多變量問題簡化為多個子問題求解，下面描述分解的子問題第k+1 次迭代求解流程。

1） 估計圖像Xk+1子問題，求解可以描述為一個典型的最小二乘模型：

通過對X一階求導(dǎo)，等價于線性方程組：

式中f ft和f ft?1分別為快速傅里葉變換和逆變換。

2） 條帶分量Sk+1子問題，條帶分量S可描述為

根據(jù)文獻[10]的加權(quán)?2,1范數(shù)最小化問題求解，定義矩陣元素的軟閾值操作符算子可求解S為

3） 變量Qk+1子問題求解，可以描述為

4） 變量Dk+1子問題求解，可以描述為

是典型的?1范數(shù)最小化問題。對?1范數(shù)最小化問題使用矩陣軟閾值操作符求解，其閉式解為

式中s hrink(l,τ)=sign(l)×max(|l|?τ,0)。

式中 κ1和κ2為收縮參數(shù)。

結(jié)合上述已解決的子問題，采用多步迭代求解，直到滿足收斂條件結(jié)束迭代循環(huán)， δ為設(shè)定的誤差參數(shù)。

2 實驗結(jié)果與分析

為了驗證WBS-MCP 模型的有效性，分別對仿真高光譜圖像和真實帶條帶的高光譜圖像去除條帶，并與UTV[7]、RBSD[10]、WLRGSD[11]和TVNLR[12]4 種性能表現(xiàn)較好的去條帶方法對比，結(jié)合視覺效果、平均峰值信噪比(MPSNR)、平均特征結(jié)構(gòu)相似性因子(MFSIM)[18]、平均等效視數(shù)(ENL)[6]和邊緣保持指數(shù)(EPI)[6]多種客觀評價方法以及譜線曲線綜合評價實驗結(jié)果。實驗測試環(huán)境為Matlab2015a，計算機硬件參數(shù)為Intel(R)Core(TM)i5-2520M CPU @2.5, RAM: 8.00 GB。

2.1 仿真數(shù)據(jù)實驗

模擬實驗中，選用了HYDICE Urban 數(shù)據(jù)集[6]和Washington DC Mall 高光譜圖像集[9]作為實驗數(shù)據(jù)源進行周期性和非周期性條帶降噪處理。為保證實驗比較的公平性，對于給定的實驗數(shù)據(jù)，每種對比方法中的參數(shù)保持一致。

2.1.1 周期性條帶實驗

圖3 多種算法對HYDICE Urban 數(shù)據(jù)集的band5 模擬仿真去條帶復(fù)原圖像視覺效果Fig.3 Restoration results of multiple algorithms using band 5 of HYDICE Urban data set

選取HYDICE Urban 數(shù)據(jù)集做實驗，添加了噪聲均值為30、比率為0.25 的周期性條帶進行性能比較。圖3 和圖4 分別為HYDICE Urban 數(shù)據(jù)集波段5 的多種算法視覺效果圖。從圖3（c）～圖3（e）中可以發(fā)現(xiàn)，在垂直方向依然存在比較明顯的條帶信息，主要原因是選取波段影像結(jié)構(gòu)突出，導(dǎo)致UTV 識別條帶區(qū)域不夠準(zhǔn)確，RBSD 在圖像變化區(qū)域的局部平滑效果不佳和WLRGSD 在影像細(xì)節(jié)區(qū)域引起影像模糊，圖3（f）的TVNLR 方法采用了非凸低秩近似，整體圖像結(jié)構(gòu)視覺效果較好，采用高階方向TV 對效果有一定提高，但由于圖像內(nèi)容豐富，無法準(zhǔn)確描述存在局部細(xì)節(jié)丟失，而本文算法采用的自適應(yīng)加權(quán)可以較好地適應(yīng)圖像內(nèi)容變化較大的情況。

圖4 多種算法對HYDICE Urban 數(shù)據(jù)集的band5 模擬仿真去條帶分量分離視覺效果Fig.4 Component results of multiple algorithms using band 5 of HYDICE Urban data set

表1 為在Washington DC Mall 高光譜圖像集中增加了5 種不同均值和比率的周期性條帶，采用不同算法去條帶MPSNR 和MSSIM 值對比。表1中粗體表示不同噪聲水平下最佳性能結(jié)果，從MPSNR 和MSSIM 指標(biāo)的比較可以看出，RBSD、UTV、WLRGSD 和TVNLR 隨著噪聲水平增加性能降低較快，特別是當(dāng)噪聲均值從25 變到30 及比率從0.25 變到0.3 時，性能指標(biāo)值降低近3 dB，而本文算法降低相對較緩慢。整體性能指標(biāo)比較結(jié)果與圖4 和圖5 的圖像質(zhì)量保持一致，特別是MSSIM指標(biāo)值在評價圖像邊緣結(jié)構(gòu)上與圖4 的視覺效果基本吻合。本文模型在多種比較方法中獲得最高的指標(biāo)性能值，這是因為本文模型從空域和光譜域聯(lián)合考慮條帶噪聲的方向特性、低秩性和結(jié)構(gòu)稀疏性等特點，特別是以波段為單位將影像分割成塊組，利用光譜的低秩特征取得更好的去條帶性能。

為了清楚地展示影像重建結(jié)果，圖5 展示了用不同方法估計Washington DC Mall 高光譜圖像集20 波段圖像的譜特征曲線。如圖5（a）所示，由于條帶和其他噪聲的影響，原始波段圖像的譜特征曲線中出現(xiàn)了許多波動?？梢杂^察到，用圖5（b）～5（f）所示的不同方法恢復(fù)后，參與比較的算法模型在一定程度上抑制了波動，圖5（b）～圖5（c）的曲線中還存在微小波動，表明圖像部分混合噪聲去除無效。從譜特征曲線結(jié)果表明，本文算法模型得到的恢復(fù)曲線明顯比其他方法平滑，說明本文算法具有更好的去噪效果。

2.1.2 非周期性條帶實驗

本實驗以Washington DC Mall 高光譜圖像部分子圖像集為數(shù)據(jù)源，加入合成的隨機分布非周期性條帶噪聲。圖6 為對加入均值為30 和比率為0.2 的隨機條帶噪聲的DC Mall 數(shù)據(jù)集波段3 的去條帶恢復(fù)圖像視覺效果。

圖6 顯示各種方法的去條帶結(jié)果和局部區(qū)域放大比較。比較圖6（c） ～圖6（e）的放大區(qū)域可以看出,恢復(fù)的影像整體偏暗，有一定的圖像模糊。從白色框圖區(qū)域可較明顯發(fā)現(xiàn)圖6（c）～圖6（f）的垂直方向結(jié)構(gòu)信息被當(dāng)作條帶去除，導(dǎo)致了不清晰,而本文算法結(jié)果基本和圖6（a）保持一致。另外，UTV 方法恢復(fù)圖像的結(jié)構(gòu)細(xì)節(jié)丟失，同時可觀察到一些殘留的條帶，與條帶方向相同的圖像結(jié)構(gòu)被抑制，缺少條帶識別能力。從圖6（d）～圖6（f）也可發(fā)現(xiàn)，由于比較方法對圖像結(jié)構(gòu)適應(yīng)性不夠，導(dǎo)致有一定的破壞，殘留少量條帶在影像中。從圖6（g）可知，本文算法恢復(fù)圖像人工痕跡不明顯，條帶被較好地去除，同時圖像結(jié)構(gòu)保留較完整。與RBSD、WLRGSD 和TVNLR 相比，本文算法利用加權(quán)的?2,1范數(shù)描述非周期條帶稀疏性，較好地保持圖像結(jié)構(gòu)，在適應(yīng)圖像結(jié)構(gòu)變化方面能力更強。

圖7 為不同去噪影像列均值曲線展示，通過列均值曲線反映去條帶噪聲影像效果，UTV 的躍變整體偏大, RBSD 和WLRGSD 在中間部分處有明顯的偏移且條帶分量偏差較多，反映出影像整體偏暗、局部區(qū)域偏亮。特別是WLRGSD 震蕩幅度過高，TVNLR 的列均值曲線也存在一定的震蕩偏差，而本文算法的列均值曲線更平滑，使得圖像整體亮度和結(jié)構(gòu)保持效果更好。

表1 多種算法對不同噪聲水平的周期條帶分量去噪MPSNR 和MFSIM 值對比Table 1 Comparison of MPSNR and MFSIM of multiple algorithms with periodic stripe component at different noise levels

圖5 對Washington DC Mall 高光譜圖像集20 波段圖像用各算法估計譜特征曲線對比Fig.5 Comparison of spectral feature curves by multiple algorithms for band 20 of Washington DC Mall

圖6 DC Mall 波段3 的非周期條帶去條帶恢復(fù)效果與局部區(qū)域放大比較Fig.6 Comparison of stripe removal restoration effect of aperiodic stripe and local area amplification of DC Mall band 3

圖7 去噪后影像列均值曲線圖Fig.7 Mean value curves of image list after denoising

2.2 真實數(shù)據(jù)實驗

為了說明本文方法的實用性，選取了兩幅受到周期噪聲和非周期噪聲污染的真實高光譜圖像集Terra MODIS進行實驗[6]，圖8 和圖9 分別給出了兩幅真實高光譜圖像集波段8 和波段34 的去條帶處理結(jié)果。

圖8 和圖9 分別為真實周期噪聲影像和真實非周期噪聲影像的去條帶視覺效果。從圖中不難發(fā)現(xiàn)，圖8（c）UTV 方法由于識別條帶能力不夠，在真實影像去除條帶時，部分條帶依然保留在影像中；RBSD 方法整體效果較好，但由于圖像自適應(yīng)加權(quán)調(diào)整不夠，依然有少量的條帶殘余沒有被去除,同時模糊了部分細(xì)節(jié)；WLRGSD 和TVNLR方法由于低秩近似忽略了圖像結(jié)構(gòu)特征，導(dǎo)致部分重要細(xì)節(jié)結(jié)構(gòu)信息被丟失，本文算法去除效果較好，主要是在考慮條帶低秩和方向特征的同時，還利用MCP 正則化約束保持結(jié)構(gòu)信息，因此

圖8 真實周期噪聲影像波段8 的去條帶視覺效果圖Fig.8 Restoration results of multiple algorithms for band 8 of real periodic noise image

本算法的去噪結(jié)果在細(xì)節(jié)和整體結(jié)構(gòu)方面保留較完整。

圖9 真實非周期噪聲影像波段34 的去條帶視覺效果圖Fig.9 Restoration results of multiple algorithms for band 34 of real aperiodic noise image

由于圖9（a）中圖像遭受不均勻細(xì)條帶干擾，從圖9 中可以看出，幾種算法去除效果有一定的差異，由于條帶分布不均勻，UTV、TVNLR和RBSD 無法適應(yīng)條帶和圖像結(jié)構(gòu)變化，在噪聲圖像中的條帶信息被大部分去除，但仍有部分殘留被當(dāng)做圖像內(nèi)容保留。WLRGSD 利用加權(quán)低秩組稀疏去除效果較好，由于本文算法采用加權(quán)?2,1范數(shù)和MCP 范數(shù)，對非周期條帶噪聲去除效果更佳。表2 為高光譜圖像客觀評價指標(biāo)對比表。

表2 高光譜圖像客觀評價指標(biāo)對比表Table 2 Comparison of objective evaluation indexes of hyperspectral images

表3 記錄了本文算法與RBSD、UTV、WLRGSD和TVNLR 4 種算法對兩個數(shù)據(jù)集中的部分波段的運行時間。從表中可以看出，由于本文算法利用了加權(quán)塊稀疏和MCP 約束，減少了部分計算量，與其他幾種算法相比，效率方面依然保持良好。

表3 多種算法運行時間對比Table 3 Comparison of running time of various algorithms s

2.3 參數(shù)分析

WBS-MCP 模型包含4 個參數(shù):圖像塊大小p、λ1、λ2、λ3。為了評估和分析這些參數(shù)的影響和最優(yōu)值，本文使用圖3 為實驗，利用PSNR 值評價測度，采用貪心策略，逐一選取參數(shù)值，選4 個參數(shù)的最優(yōu)值。從表4 中不同圖像塊取值下的PSNR 值對比可以發(fā)現(xiàn)，當(dāng)p取值為7 時，獲得最佳的性能，隨著圖像塊加大，性能有一定下降。因此，本文選取p=7。

表4 不同p 值的PSNR 指標(biāo)值對比Table 4 PSNR comparison of different values of p

圖10 繪制了不同正則化參數(shù)下λ1、λ2、λ3的PSNR 值實驗結(jié)果。圖10（a）隨著λ1從0.001 增加到0.005，PSNR 值有明顯的提高，λ1最佳取值在0.005 附近。同理，圖10（b）中λ2最佳取值在0.045 附近, 圖10（c）中λ3最佳取值在0.002 附近。

圖10 不同正則化參數(shù)下λ1、λ2、λ3 的PSNR 曲線圖Fig.10 PSNR curves of λ1， λ2， λ3 with different regularization parameters

3 結(jié)論

本文提出一種針對高光譜圖像直條帶去除的新模型WBS-MCP。在分析了條帶噪聲的特點基礎(chǔ)上，該模型利用MCP 范數(shù)和加權(quán)?2，1范數(shù)約束，較好地刻畫條帶低秩特性，同時具有經(jīng)典UTV正則化項模型，能準(zhǔn)確提取圖像的主要結(jié)構(gòu)信息，采用自適應(yīng)的權(quán)值優(yōu)化了圖像結(jié)構(gòu)約束，較好地保留估計影像的邊緣和紋理結(jié)構(gòu)等信息。通過在模擬數(shù)據(jù)和真實數(shù)據(jù)的周期性和非周期條帶去除，并采用多種評價指標(biāo)與同類性能優(yōu)秀的算法比較，可以看出，WBS-MCP 模型采用了非凸MCP 來表征條帶噪聲的低秩性，可以有效地將條帶分量和圖像分離，同時加權(quán)?2，1范數(shù)約束在圖像結(jié)構(gòu)保持方面都優(yōu)于其他比較方法。將來的工作和研究方向是嘗試將旋轉(zhuǎn)策略加入到WBS-MCP 模型中來捕獲斜條帶。此外，還可以將本方法擴展到視頻雨斑去除和紅外圖像條帶噪聲去除等應(yīng)用中。