CN106097278B - 一种多维信号的稀疏模型、重建方法和字典训练方法 - Google Patents

一种多维信号的稀疏模型、重建方法和字典训练方法 Download PDF

Info

Publication number
CN106097278B
CN106097278B CN201610474921.0A CN201610474921A CN106097278B CN 106097278 B CN106097278 B CN 106097278B CN 201610474921 A CN201610474921 A CN 201610474921A CN 106097278 B CN106097278 B CN 106097278B
Authority
CN
China
Prior art keywords
sparse
dictionary
tensor
training
dimensional
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201610474921.0A
Other languages
English (en)
Other versions
CN106097278A (zh
Inventor
齐娜
施云惠
尹宝才
丁文鹏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing University of Technology
Original Assignee
Beijing University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing University of Technology filed Critical Beijing University of Technology
Priority to CN201610474921.0A priority Critical patent/CN106097278B/zh
Publication of CN106097278A publication Critical patent/CN106097278A/zh
Application granted granted Critical
Publication of CN106097278B publication Critical patent/CN106097278B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/70—Denoising; Smoothing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20081—Training; Learning

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Image Analysis (AREA)

Abstract

本发明公开一种多维信号的稀疏模型，其能够保证不需要采用Kronecker乘积，从而保证在算法复杂度和存储空间上都有明显的改进。这种多维信号的稀疏模型，其为公式
Figure DDA0001959402860000011
其中，张量X表示成一个N维稀疏张量
Figure DDA0001959402860000012
与一系列稀疏字典
Figure DDA0001959402860000013
的张量乘积，In≤Mn，Dn定义为第n个方向的字典，K是稀疏度，用来刻画稀疏系数

Description

一种多维信号的稀疏模型、重建方法和字典训练方法
技术领域
本发明属于稀疏表示的稀疏重建和字典训练的技术领域，具体地涉及一种多维信号的稀疏模型、重建方法和字典训练方法。
背景技术
稀疏表示作为近年来图像视频建模的有效方法，成功地应用于诸如图像去噪、超分辨率重建、人脸识别等计算机视觉领域。传统的稀疏表示通常是将多维信号转化为一维信号，该一维信号可以由一个字典中的几个基元线性组合表示。高阶信号(图像、视频等信号)首先需要被转化为一维信号，利用一些处理向量的方法对其处理。传统稀疏表示模型的研究包括模型的建立，稀疏重建，和字典训练方法。尽管一维稀疏表示模型有效，但是由于在处理高维信号时，随之维度的增加，向量化后的一维信号是一个高维向量，在其进行处理时，需要较大维度的观测矩阵，会引起维度灾难等问题，因此存在空间存储和算法复杂度的压力。
目前的一些研究表明，在对多维信号进行处理时，保持多维信号本身的空间结构存在很多好处。针对二维图像信号，提出的二维稀疏表示模型充分利用图像内部的二维空间结构和相关性，被成功地应用于图像去噪和图像超分辨率重建。其中二维字典训练问题可以通过两阶段的块坐标松弛方法来实现。二维稀疏重建问题可以通过诸如2DSL0或2D-OMP等一维稀疏重建的拓展方法来解决，也可以利用克罗奈克积(Kronecker product)将其转化为一维稀疏重建问题，通过一些常用的一维稀疏重建方法求解。二维稀疏表示模型的采用可以大幅度地减少图像处理的时间复杂度和存储空间。但是由于需要利用一维稀疏重建方法，将二维稀疏表示模型直接推广到多维稀疏表示模型是比较困难的。
此外，张量结构也被引入到向量的稀疏表示中，用来逼近多维信号在每一个模式(方向)的结构。由于有约束的Tucker模型和张量的Kronecker表达具有等价性，同时考虑到多维稀疏特性和块稀疏特性等稀疏约束，张量可以由给定的不同方向的可分离Kronecker字典表示。相应的Kronecker OMP和N-way Block OMP方法用来解决字典给定情况下的稀疏重建问题，由于需要使用Kronecker乘积运算，这两种方法都相对比较复杂。另外，还有一类利用张量结构的模型是基于张量分解的，如PARAFAC，基于Tucker模型的张量分解；或者是张量低秩逼近，如LRTA，HOSVD，和Tensor-DL等。然而这些基于张量的方法都是对张量信号本身进行分解或逼近，目前还没有统一的框架来刻画张量的稀疏表示模型，在算法复杂度和存储空间上存在矛盾。
发明内容
本发明的技术解决问题是：克服现有技术的不足，提供一种多维信号的稀疏模型，其能够保证不需要采用Kronecker乘积，从而保证在算法复杂度和存储空间上都有明显的改进。
本发明的技术解决方案是：这种多维信号的稀疏模型，其为公式(1)
Figure GDA0001959402850000021
其中，张量
Figure GDA0001959402850000022
表示成一个N维稀疏张量
Figure GDA0001959402850000023
与一系列稀疏字典
Figure GDA0001959402850000024
的张量乘积，In≤Mn，Dn定义为第n个方向的字典，K是稀疏度，用来刻画稀疏系数
Figure GDA0001959402850000025
中的非零元素的个数。
还提供了一种多维信号的稀疏模型的重建方法，该方法为基于张量的迭代收缩阈值方法TISTA，针对公式(3)和(2)得到公式(10)、(11)
Figure GDA0001959402850000031
Figure GDA0001959402850000032
其中，Sτ(·)是软阈值算子，定义为sign(·)max(|·|-τ,0)，而相应地，Hτ(·)是硬阈值算子，定义为max(|·|-τ,0)；
Figure GDA0001959402850000033
ηk≥1；
Figure GDA0001959402850000034
满足公式(7)
Figure GDA0001959402850000035
不论是约束是
Figure GDA0001959402850000036
还是
Figure GDA0001959402850000037
将上述算子记做公式(12)
Figure GDA0001959402850000038
此外，通过考虑
Figure GDA0001959402850000039
来加速上述迭代收敛过程，其中ξk>0是一个合适的步长，
Figure GDA00019594028500000310
即实际求解公式
Figure GDA00019594028500000311
还提供了重建方法的字典训练方法，字典训练模型为公式(13)
Figure GDA00019594028500000312
其中
Figure GDA00019594028500000313
是张量样本集I中所有训练样本的稀疏系数集合，λ是用来平衡保真度和稀疏度的，其中
Figure GDA00019594028500000314
为待训练的可分离的字典。
本发明旨在利用张量对多维信号建立一种新的稀疏表示模型，即张量稀疏表示模型，同时提出其相应的多维稀疏重建和多维字典训练方法。提出的张量稀疏表示模型能够通过字典学习得到一系列自适应的可分离字典，来更好地刻画多维信号内部每个方向的结构。本发明提出的基于最近邻方法的多维稀疏重建算法能够更进一步地降低算法复杂度。本发明通过三维多光谱图像的去噪验证了提出的张量稀疏表示模型能够大幅度地降低算法复杂度和存储空间，但是却可以获得和最先进的稀疏表示模型相当的去噪性能。因此，张量稀疏表示模型是一种更好地刻画多维信号的稀疏表达能力的方式。
附图说明
图1示出了根据本发明的基于张量的迭代收缩阈值方法的流程图；
图2示出了根据本发明的基于张量的多维字典训练方法的流程图。
图3示出了根据本发明的多光谱图像去噪的流程图。
具体实施方式
这种多维信号的稀疏模型，其为公式(1)
Figure GDA0001959402850000041
其中，张量
Figure GDA0001959402850000042
表示成一个N维稀疏张量
Figure GDA0001959402850000043
与一系列稀疏字典
Figure GDA0001959402850000044
的张量乘积，In≤Mn，Dn定义为第n个方向的字典，K是稀疏度，用来刻画稀疏系数
Figure GDA0001959402850000045
中的非零元素的个数。
优选地，给定稀疏字典
Figure GDA0001959402850000046
相应的稀疏模型为公式(2)
Figure GDA0001959402850000047
其中，λ是用来平衡保真度和稀疏度的。
优选地，通过松弛转化为l1约束的凸规划问题，相应的稀疏模型为公式(3)
Figure GDA0001959402850000048
还提供了一种多维信号的稀疏模型的重建方法，该方法为基于张量的迭代收缩阈值方法TISTA(如图1所示)，针对公式(3)和(2)得到公式(10)、(11)
Figure GDA00019594028500000515
Figure GDA0001959402850000051
其中，Sτ(·)是软阈值算子，定义为sign(·)max(|·|-τ,0)，而相应地，Hτ(·)是硬阈值算子，定义为max(|·|-τ,0)；
Figure GDA0001959402850000052
ηk≥1；
Figure GDA0001959402850000053
满足公式(7)
Figure GDA0001959402850000054
不论是约束是
Figure GDA0001959402850000055
还是
Figure GDA0001959402850000056
将上述算子记做公式(12)
Figure GDA0001959402850000057
优选地，通过考虑
Figure GDA0001959402850000058
来加速上述迭代收敛过程，即实质上是采用如下公式
Figure GDA0001959402850000059
其中ξk>0是一个合适的步长，
Figure GDA00019594028500000510
Figure GDA00019594028500000511
还提供了重建方法的字典训练方法(如图2所示)，字典训练模型为公式(13)
Figure GDA00019594028500000512
其中
Figure GDA00019594028500000513
是张量样本集I中所有训练样本的稀疏系数集合，λ是用来平衡保真度和稀疏度的，其中
Figure GDA00019594028500000514
为待训练的可分离的字典。
优选地，稀疏重建问题是解决在给定字典
Figure GDA0001959402850000061
下，训练样本集的
Figure GDA0001959402850000062
的稀疏系数
Figure GDA0001959402850000063
其目标函数为
Figure GDA0001959402850000064
该问题通过基于张量的迭代收缩阈值方法TISTA进行求解。
优选地，当更新字典Dn时，通过
Figure GDA0001959402850000065
计算得到
Figure GDA0001959402850000066
以实现
Figure GDA0001959402850000067
然后按第n模式方式展开
Figure GDA0001959402850000068
得到
Figure GDA0001959402850000069
从而使得
Figure GDA00019594028500000610
最后Dn通过解决公式(15)的二次约束二次规划问题来求解：
Figure GDA00019594028500000611
其中
Figure GDA00019594028500000612
分别为
Figure GDA00019594028500000613
按n-模式方向展开的矩阵，该问题通过Lagrange对偶方法求解，最终
Figure GDA00019594028500000614
其中Λ是对偶变量组成的对偶矩阵。
本发明旨在利用张量对多维信号建立一种新的稀疏表示模型，即张量稀疏表示模型，同时提出其相应的多维稀疏重建和多维字典训练方法。提出的张量稀疏表示模型能够通过字典学习得到一系列自适应的可分离字典，来更好地刻画多维信号内部每个方向的结构。本发明提出的基于最近邻方法的多维稀疏重建算法能够更进一步地降低算法复杂度。本发明通过三维多光谱图像的去噪验证了提出的张量稀疏表示模型能够大幅度地降低算法复杂度和存储空间，但是却可以获得和最先进的稀疏表示模型相当的去噪性能。因此，张量稀疏表示模型是一种更好地刻画多维信号的稀疏表达能力的方式。
下面更详细地说明本方法的实施例。
传统的稀疏表示模型是基于向量形式的数据进行建模和表达的。它将一维信号向量化，得到的向量形式的一维信号可以表达为一个大的字典下少量基元的线性组合。一维稀疏表示模型尽管简单，但是它忽略了多维信号的内部结构和破坏了其相关性。将多维信号转化为一维信号会生成高维信号，对于处理高维信号而言，会大大增加存储需求和计算资源。因此，本发明旨在利用张量对多维信号建立一种新的稀疏表示模型，即张量稀疏表示模型，同时提出其相应的多维稀疏重建和多维字典训练方法。提出的张量稀疏表示模型能够通过字典学习得到一系列自适应的可分离字典，来更好地刻画多维信号内部每个方向的结构。本发明提出的基于最近邻方法的多维稀疏重建算法能够更进一步地降低算法复杂度。本发明通过三维多光谱图像的去噪验证了提出的张量稀疏表示模型能够大幅度地降低算法复杂度和存储空间，但是却可以获得和最先进的稀疏表示模型相当的去噪效果。因此，张量稀疏表示模型是一种更好地刻画多维信号的稀疏表达能力的方式。
为了便于下文中公式和符号的方便理解，这里先给出些符号公式的解释。N阶张量记作
Figure GDA0001959402850000071
N阶张量
Figure GDA0001959402850000072
的l0，l1，lF范数分别定义如下：
Figure GDA0001959402850000073
Figure GDA0001959402850000074
其中
Figure GDA0001959402850000075
是张量
Figure GDA0001959402850000076
的第
Figure GDA0001959402850000077
个元素。通过固定其他方向的索引，而遍历某一个方向(第n个方向)的索引，可以得到按第n模式方向展开的n模式向量。张量
Figure GDA0001959402850000078
按第n模式展开得到矩阵
Figure GDA0001959402850000079
定义为将所有的n模式向量按列排列为一个矩阵。张量
Figure GDA00019594028500000710
与矩阵
Figure GDA00019594028500000711
的按n模式的乘积定义为
Figure GDA00019594028500000712
可以得到一个N阶张量
Figure GDA00019594028500000713
它的每一个元素为
Figure GDA00019594028500000714
可以通过公式
Figure GDA00019594028500000715
得到。两个相同大小的张量的内积定义为它们中每一个元素对应位置相乘的和，记作
Figure GDA00019594028500000716
符号
Figure GDA00019594028500000717
表示克洛内克积(Kronecker乘积)。对张量
Figure GDA00019594028500000718
向量化可以得到向量x。
具体地，给定一个N维信号表示为张量
Figure GDA00019594028500000719
为了能够刻画其结构以及表达其不同维度的相关性，提出基于张量的稀疏重建模型如下
Figure GDA0001959402850000081
这里，张量
Figure GDA0001959402850000082
表示成一个N维稀疏张量
Figure GDA0001959402850000083
与一系列稀疏字典
Figure GDA0001959402850000084
的张量乘积，In≤Mn，Dn定义为第n个方向的字典，K是稀疏度，用来刻画稀疏系数
Figure GDA0001959402850000085
中的非零元素的个数。
针对张量的稀疏表示模型，以下从两方面稀疏重建问题和字典训练问题的来介绍本发明所提出的求解方法。
给定稀疏字典
Figure GDA0001959402850000086
相应的稀疏重建模型为，
Figure GDA0001959402850000087
其中，λ是用来平衡保真度和稀疏度的。此外，这个l0约束的非凸优化问题，可以通过松弛转化为l1约束的凸规划问题，定义如下：
Figure GDA0001959402850000088
本发明提出了一种基于张量的迭代收缩阈值方法(Tensor-based IterativeShrinkage Thresholding Algorithm,TISTA)的方法来解决上述稀疏重建问题，主要是通过迭代收缩方法来实现。
具体如下，上述问题可记作：
Figure GDA0001959402850000089
其中
Figure GDA00019594028500000810
表示数据拟合项
Figure GDA00019594028500000811
而
Figure GDA00019594028500000812
表示稀疏约束项，即
Figure GDA00019594028500000813
或
Figure GDA00019594028500000814
上述非光滑的最优化问题可以通过迭代收缩方法来进行求解，
Figure GDA00019594028500000815
写做是f在附近的
Figure GDA00019594028500000816
进行线性化逼近，同时考虑最近邻正则项和非光滑正则项，则
Figure GDA00019594028500000817
可以通过如下最优化问题来更新：
Figure GDA00019594028500000818
这里Lk>0是一个Lipschitz常数，而
Figure GDA00019594028500000819
是定义在张量场的梯度。上述问题可以转化为：
Figure GDA00019594028500000820
导出
Figure GDA00019594028500000821
在张量场的梯度
Figure GDA00019594028500000822
满足
Figure GDA0001959402850000091
接下来我们讨论Lk>0的求解。假设f是一个C1,1的光滑凸函数。对于任意的
Figure GDA0001959402850000092
f是连续可导的，且具有Lipschitz常数的连续梯度，L(f)满足，
Figure GDA0001959402850000093
这里的||·||表示定义在张量上的F范数，L(f)是f的Lipschitz常数，通过推导，可以得到，
Figure GDA0001959402850000094
因此，梯度
Figure GDA0001959402850000095
的最小Lipschitz常数是
Figure GDA0001959402850000096
因此在具体过程中，
Figure GDA0001959402850000097
最后，在给定不同的稀疏度约束
Figure GDA0001959402850000098
下，本发明给出了不同的稀疏重建解法。针对问题(3)和(2)，可以得到如下解
Figure GDA0001959402850000099
Figure GDA00019594028500000910
其中，Sτ(·)是软阈值算子，定义为sign(·)max(|·|-τ,0)，而相应地，Hτ(·)是硬阈值算子，定义为max(|·|-τ,0)。最后，不论是约束是
Figure GDA00019594028500000911
还是
Figure GDA00019594028500000912
可以将上述算子记做
Figure GDA00019594028500000913
此外，可以通过考虑
Figure GDA00019594028500000914
而不是
Figure GDA00019594028500000915
来加速上述迭代收敛过程，其中ξk>0是一个合适的步长。
Figure GDA00019594028500000916
这里的
Figure GDA00019594028500000917
这种外插的方式类似于在FISTA算法中的加速方法。表1总结了本发明提出来的TISTA算法。
表1
Figure GDA0001959402850000101
下面介绍本发明提出来的多维字典训练问题的方法。
相应地，对于给定一个张量的训练样本集合
Figure GDA0001959402850000102
包括S个张量
Figure GDA0001959402850000103
本发明提出的基于张量的字典训练模型如下：
Figure GDA0001959402850000104
其中
Figure GDA0001959402850000105
是张量样本集I中所有训练样本的稀疏系数集合，λ是用来平衡保真度和稀疏度的，其中
Figure GDA0001959402850000106
为待训练的可分离的字典。
该问题可以通过两阶段——稀疏重建和字典更新的分块松弛方法来迭代求解，直到达到预设的迭代终止条件为止才停止循环，例如目标函数的相对误差小于预设的误差。表2总结了基于张量的字典训练方法：
具体地，稀疏重建问题是解决在给定字典
Figure GDA0001959402850000107
下，训练样本集的
Figure GDA0001959402850000111
的稀疏系数
Figure GDA0001959402850000112
其目标函数为
Figure GDA0001959402850000113
该问题可以直接通过上述本发明提出的TISTA的稀疏重建算法进行求解，而不是通过Kronecker product将其转化为一维稀疏重建问题来求解。此外，我们这里的目标函数更加灵活，可以同时求解整个训练样本集I的稀疏系数，也可以先将I分段，然后分段求解训练样本集中的稀疏系数，也可以并行计算其稀疏系数。因此，这种稀疏重建算法更加地灵活。
表2
Figure GDA0001959402850000114
字典更新过程主要是在利用求解得到的稀疏系数来更新字典
Figure GDA0001959402850000115
更新字典
Figure GDA0001959402850000116
的最优化过程相类似，不失一般性，以更新字典Dn为例来介绍本发明的字典训练方法。
由于在张量稀疏模型中n-mode积(n模式乘积)的可交换性，基于本发明提出的张量稀疏表示模型，每一个样本满足
Figure GDA0001959402850000117
其中
Figure GDA0001959402850000121
因此很容易通过按第n模式方向展开
Figure GDA0001959402850000122
得到
Figure GDA0001959402850000123
而不是通过利用克罗内克积Kronecker积的方式实现(根据张量运算，第n-mode展开的矩阵满足
Figure GDA0001959402850000124
Figure GDA0001959402850000125
因此，当更新字典Dn时，首先通过
Figure GDA0001959402850000126
计算得到
Figure GDA0001959402850000127
以实现
Figure GDA0001959402850000128
然后按第n模式方式展开
Figure GDA0001959402850000129
得到
Figure GDA00019594028500001210
从而使得
Figure GDA00019594028500001211
最后Dn可以通过解决如下的二次约束二次规划问题来求解：
Figure GDA00019594028500001212
其中
Figure GDA00019594028500001213
分别为
Figure GDA00019594028500001214
按n-模式方向展开的矩阵。该问题可以通过Lagrange对偶方法求解，最终
Figure GDA00019594028500001215
其中Λ是对偶变量组成的对偶矩阵。这种求解
Figure GDA00019594028500001216
避免了使用Kronecker积，可以更加有效地降低更新字典的算法复杂度。
为了说明本发明的有效性，本发明将提出的张量稀疏模型，以及其相应的稀疏重建方法和字典训练方法应用于多光谱图像去噪问题。
具体地，利用哥伦比亚大学的多光谱图像数据集进行测试。该图像数据集有32个真实场景，每一个场景的空间分辨率是512×512，频谱为31个谱带，谱带范围是从400nm到700nm。对多光谱图像添加不同的均方误差(σ＝5,10,20,30,50)的高斯噪声后得到噪声图像，然后对噪声图像做去噪处理。本发明提出的基于张量稀疏模型的方法中，三个方向的字典D1,D2,D3的大小都是5×10，同时利用冗余离散余弦字典(Overcomplete DictionaryLearning,ODCT)进行字典的初始化，针对某一副噪声图像，从该图像中随机采样40,000块大小为5×5×5的立方体(cubes)，然后通过本发明提出来的多维字典训练方法训练得到字典，通过训练得到的字典，对噪声多光谱图像采样得到的所有立方体块，利用本发明提出的多维稀疏重建方法进行稀疏系数的求解，从而重建去噪的MSI图像块，最终实现多光谱图像的去噪。其中，对应不同的均方误差σ＝5,10,20,30,50，相应地控制稀疏度约束的参数λ＝9,20,45,70,160。
为了验证以上多光谱图像的每一个谱带的图像的去噪效果，主要是通过峰值信噪比(Peak Signal to Noise Ratio，PSNR)度量，单位为分贝(dB)。
其计算公式如下：
Figure GDA0001959402850000131
两幅大小为m*n的图像的均方误差MSE的定义如下：
Figure GDA0001959402850000132
其中I,J分别表示原始不含噪声的图像和利用稀疏编码方法重建图像，而I(x,y),J(x,y)为对应于位置(x,y)处的像素值，则均方误差越小，则PSNR越高，则该方法的去噪效果越高。
客观质量的评判标准为结构相似度量：它是一种基于结构失真图像质量的评价方法，该方法将亮度和对比度从图像结构信息中分离，并结合结构信息对图像质量进行评价.定义了SSIM：
SSIM(i,j)＝[L(i,j)]α·[C(i,j)]β·[S(i,j)]γ (18)
其中：
Figure GDA0001959402850000133
Figure GDA0001959402850000134
Figure GDA0001959402850000135
其中L表示亮度(Lightness),用均值(μi,μj)作为亮度的估计，C表示对比度(Contrast),用标准差(σi,σj)作为对比度的估计，而相应的协方差σij作为结构相似程度的度量。而α,β,γ是用来调整亮度，对比度和结构信息的权重，为了防止分母出现零或接近0，而产生不稳定现象，而引入了λ1,λ2,λ3。当α＝β＝γ＝1，λ3＝λ2/2时，公式简化为
Figure GDA0001959402850000136
由于SSIM的较好的特性：对称性，有界性，唯一性。使得它很好的反映了图像的客观质量。尤其是由于观察值在某一瞬间对某个局部区域的细节更关注些，所以图像的整体统计特征不能反映图像的特点。因此SSIM结合PSNR更能反映图像的质量。一般PSNR较大时SSIM的区分度较小，即PSNR高，SSIM也高，而在PSNR较小时，SSIM就具有较好的区分度。因此这里主要是通过PSNR和SSIM来评价超分辨图像的结果。
表3给出了多光谱图像在不同的均方误差情况下的平均峰值信噪比(PSNR)和结构相似度量(SSIM)。本发明主要和4个方法进行了比较，包括PARAFAC,Bw-KSVD，3DK-SVD和基于低秩张量逼近的方法LRTA。由表3可以看出本文的方法优于其他三种方法(PARAFAC,Bw-KSVD，3DK-SVD),与LRTA的结果相当。
表3
Figure GDA0001959402850000141
这里需要指出的是，在进行上述实验过程中，与基于稀疏表示的方法(Bw-KSVD,3DK-SVD)相比，本发明的字典大小是最小的。本发明的字典大小为3×5×10，其中每个方向的字典大小为5×10，而Bw-KSVD方法的字典大小分别为64×128×31，其中31副多光谱图像的每一谱带都需要一个64×128大小的字典。3DK-SVD方法所需要的字典大小是448×500，其中每一个测试的基本图像块单位为一个8×8×7的大小。因此Bw-KSVD方法和3DK-SVD方法所需字典是本发明方法字典的1693，1493倍。
此外，为了说明本发明的算法复杂度低，我们同时给出了本发明中TISTA算法和传统的一维稀疏重建方法——FISTA算法的复杂度比较。具体地，我们从一些局部的多光谱图像(Cube Size大小分别为：12×12×31，16×16×31，32×32×31)中采样一些5×5×5的立体块Cubes，Cubes的个数分别为1758，3888，21168。然后设置三个方向的字典D1,D2,D3的大小都是5×10，且都为冗余离散余弦字典(Overcomplete Dictionary Learning,ODCT)。同时设置迭代次数num＝50，λ＝1。由于TISTA算法可以单独考虑每一个样本的重建(single)，也可以成批地考虑(Batch)或者是所有的样本(All)一起考虑。表4所示给出了求解每个集合的稀疏重建问题在不同情况下所需要的时间(秒,seconds)，另外也给出了传统的FISTA方法的稀疏重建问题所需的时间。通过表4可以看出本发明提出的TISTA方法比FISTA所用的时间要少很多。因此本发明的方法在算法复杂度上有明显地提升。
表4
Figure GDA0001959402850000151
这里首先介绍给在给定一个张量信号
Figure GDA0001959402850000152
和一系列的各个方向的字典
Figure GDA0001959402850000153
时，如何求解张量信号
Figure GDA0001959402850000154
的在给定控制稀疏度约束的参数λ的条件下的稀疏系数的问题。
1.输入一个张量信号
Figure GDA0001959402850000155
和一系列的各个方向的字典
Figure GDA0001959402850000156
2.初始化，设置待估计的稀疏系数为
Figure GDA0001959402850000157
迭代步长t1＝1。本发明中这里是随机初始化稀疏系数
Figure GDA0001959402850000158
3.计算Lipschitz常数
Figure GDA0001959402850000159
这里我们将ηk设置为1。
4.计算f在最近邻点
Figure GDA00019594028500001510
的梯度。
Figure GDA00019594028500001511
5.利用Proximal算子求解稀疏系数
Figure GDA00019594028500001512
即
Figure GDA00019594028500001513
Figure GDA00019594028500001514
Figure GDA00019594028500001515
其中，Sτ(·)是软阈值算子，定义为sign(·)max(|·|-τ,0)，而相应地，Hτ(·)是硬阈值算子，定义为max(|·|-τ,0)。
6.计算迭代步长
Figure GDA00019594028500001516
其中
Figure GDA00019594028500001518
7.更新最近邻稀疏系数
Figure GDA00019594028500001519
8.判断是否达到迭代次数，未达到迭代次数回到第三步，已达到则直接输出当前的稀疏系数
Figure GDA0001959402850000161
再介绍当得到一组张量信号组成的训练样本集后，在给定稀疏约束的条件下，训练得到该训练样本的多维稀疏字典
Figure GDA0001959402850000162
具体如下：
1.构造训练样本集
Figure GDA0001959402850000163
其中S为样本个数，每个样本
Figure GDA0001959402850000164
大小为
Figure GDA0001959402850000165
在这个具体的多光谱图像中，具体的样本
Figure GDA0001959402850000166
大小为5×5×5的立方体，从多光谱图像得到的所有的立方体中，随机选择40,000，作为训练样本。
2.初始化各个方向的字典
Figure GDA0001959402850000167
针对多光谱图像，我们的字典为D1,D2,D3，采用冗余离散余弦字典(Overcomplete Dictionary Learning,ODCT)进行字典的初始化。
3.多维稀疏求解阶段，根据当前给定的多维字典，求解所有样本的稀疏系数。具体地，利用本发明提出的基于张量的迭代收缩阈值方法求解如下稀疏重建问题。
Figure GDA0001959402850000168
4.多维字典更新阶段。依次更新字典
Figure GDA0001959402850000169
具体地，针对Dn的更新，首先计算
Figure GDA00019594028500001610
然后将
Figure GDA00019594028500001611
按第n个模式展开得到相应的矩阵
Figure GDA00019594028500001612
则具体地，可以通过如下方程更新字典Dn。
Figure GDA00019594028500001613
其他方向的字典
Figure GDA00019594028500001614
的更新过程类似。
5.判断是否达到迭代终止条件，如1)迭代条件为迭代次数是否达到上限num，2)迭代条件为噪声误差是否达到指定值。
6.输出训练得到的字典
Figure GDA00019594028500001615
基于张量稀疏表示模型的多维信号去噪的实施例
1.利用含噪的多光谱图像构造待求解的多维信号，每个多维信号大小为5×5×5，共采样N块，得到待重建的多光谱图像块的集合In。
2.利用上述基于张量稀疏表示模型的多维字典训练方法训练得到的字典
Figure GDA00019594028500001616
针对待重建的多光谱图像块的集合，利用基于张量的迭代收缩阈值方法求解稀疏系数。
3.根据稀疏系数和字典
Figure GDA0001959402850000171
可以重建多光谱图像块，得到重建图像块集合Ire。
4.利用上述N个多光谱图像块重建去噪多光谱图像。
根据采样过程中的采样方式以及相应图像块的重叠方式，将目前得到的N个重建图像块重新恢复回原始图像大小，对应于有重叠的地方，则采用取平均值的操作。
以上所述，仅是本发明的较佳实施例，并非对本发明作任何形式上的限制，凡是依据本发明的技术实质对以上实施例所作的任何简单修改、等同变化与修饰，均仍属本发明技术方案的保护范围。

Claims (3)

1.一种多光谱图像去噪方法，其特征在于包括以下步骤：
(a)利用含噪的多光谱图像构造待求解的多维信号，每个多维信号大小为I1×I2×…×In…×IN，共采样O块，得到待重建的多光谱图像块的集合In；
(b)利用基于张量稀疏表示模型的多维字典训练方法训练得到的字典
Figure FDA0003277004740000011
针对待重建的多光谱图像块的集合，利用基于张量的迭代收缩阈值方法求解稀疏系数；
(c)根据稀疏系数和训练得到的字典
Figure FDA0003277004740000012
可以重建多光谱图像块，得到重建图像块集合Ire；
(d)利用步骤(c)的O个重建图像块重建去噪多光谱图像；
(e)根据采样过程中的采样方式以及相应图像块的重叠方式，将目前得到的O个重建图像块重新恢复回原始图像大小，对于有重叠的地方，则采用取平均值的操作；
步骤(b)中，字典训练模型为公式(13)
Figure FDA0003277004740000013
Figure FDA0003277004740000014
其中
Figure FDA0003277004740000015
是张量样本集I中所有训练样本的稀疏系数集合，
Figure FDA0003277004740000016
为训练样本集I中
Figure FDA0003277004740000017
对应的稀疏系数，
Figure FDA0003277004740000018
的定义是实数域
Figure FDA0003277004740000019
上大小为M1×M2×…×Mn…×MN×S的N+1阶张量；
Figure FDA00032770047400000110
为待训练的可分离的字典,其中
Figure FDA00032770047400000111
Mn为稀疏字典矩阵Dn的列数，In为稀疏字典矩阵Dn的行数，r是稀疏字典Dn的列的索引；×n表示张量按n模式的乘积；λ是用来平衡保真度和稀疏度的系数；字典训练模型公式(13)的求解通过解决稀疏重建问题和字典训练问题交替迭代来实现的，包括以下步骤：
(I)构造训练样本集
Figure FDA0003277004740000021
其中S为样本个数，每个样本
Figure FDA0003277004740000022
大小为
Figure FDA0003277004740000023
(II)初始化各个方向的字典
Figure FDA0003277004740000024
针对多光谱图像，采用冗余离散余弦字典进行字典的初始化；
(III)多维稀疏重建阶段，根据当前给定的多维字典，求解所有样本的稀疏系数，利用基于张量的迭代收缩阈值方法求解稀疏重建问题；
(IV)多维字典更新阶段，依次更新字典
Figure FDA0003277004740000025
(V)判断是否达到迭代终止条件；未达到迭代条件回到步骤(III)，已达到则继续下一步；
(VI)输出训练得到的字典
Figure FDA0003277004740000026
2.根据权利要求1所述的多光谱图像去噪方法，其特征在于所述基于张量的迭代收缩阈值方法，包括以下步骤：
(1)输入一个张量信号χ，以及一系列的各个方向的字典
Figure FDA0003277004740000027
(2)初始化，设置待估计的稀疏系数为
Figure FDA0003277004740000028
迭代步长t1＝1，其中
Figure FDA0003277004740000029
为初始化的稀疏系数；
(3)计算Lipschitz常数
Figure FDA00032770047400000210
ηk设置为1；
(4)根据公式(7)计算f在最近邻点Ck-1的梯度
Figure FDA00032770047400000211
其中，f的定义是：输入一个张量信号
Figure FDA00032770047400000212
以及一系列的各个方向的字典
Figure FDA00032770047400000213
求解其稀疏系数
Figure FDA00032770047400000214
Figure FDA00032770047400000215
表示数据拟合项
Figure FDA0003277004740000031
Figure FDA0003277004740000032
写作是f在附近的
Figure FDA0003277004740000033
进行线性化逼近，可以表达为
Figure FDA0003277004740000034
其中
Figure FDA0003277004740000035
为利用基于张量的迭代收缩阈值方法TISTA求解
Figure FDA0003277004740000036
的稀疏系数时第k-1次得到的稀疏系数，k表示迭代次数；
(5)利用公式(12)求解稀疏系数
Figure FDA0003277004740000037
Figure FDA0003277004740000038
其中，
Figure FDA0003277004740000039
是阈值算子，当稀疏约束为l1范数约束的话，阈值算子
Figure FDA00032770047400000310
为软阈值；软阈值算子Sτ(·)定义为sign(·)max(|·|-τ,0)，其中τ为阈值参数，
Figure FDA00032770047400000311
则公式(12)表达式为：
Figure FDA00032770047400000312
当稀疏约束为l0范数约束的话，阈值算子
Figure FDA00032770047400000313
为硬阈值；硬阈值算子Hτ(·)定义为max(|·|-τ,0)；则公式(12)表达式为：
Figure FDA00032770047400000314
(6)计算迭代步长
Figure FDA00032770047400000315
其中
Figure FDA00032770047400000316
tk是第k次的步长；
(7)更新最近邻稀疏系数
Figure FDA00032770047400000317
其中
Figure FDA00032770047400000318
为求解
Figure FDA00032770047400000319
的稀疏系数时第k-2次得到的稀疏系数；
(8)判断是否达到迭代次数，未达到迭代次数回到步骤(3)，已达到则直接输出当前的稀疏系数
Figure FDA00032770047400000320
3.根据权利要求2所述的多光谱图像去噪方法，其特征在于：所述步骤(IV)中，针对Dn的更新，首先计算
Figure FDA0003277004740000041
然后将
Figure FDA0003277004740000042
按第n个模式展开得到相应的矩阵
Figure FDA0003277004740000043
最后Dn通过解决公式(15)的二次约束二次规划问题来求解：
Figure FDA0003277004740000044
其中
Figure FDA0003277004740000045
分别为
Figure FDA0003277004740000046
按n-模式方向展开的矩阵，该问题通过Lagrange对偶方法求解，最终
Figure FDA0003277004740000047
其中Λ是对偶变量组成的对偶矩阵。
CN201610474921.0A 2016-06-24 2016-06-24 一种多维信号的稀疏模型、重建方法和字典训练方法 Active CN106097278B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610474921.0A CN106097278B (zh) 2016-06-24 2016-06-24 一种多维信号的稀疏模型、重建方法和字典训练方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610474921.0A CN106097278B (zh) 2016-06-24 2016-06-24 一种多维信号的稀疏模型、重建方法和字典训练方法

Publications (2)

Publication Number Publication Date
CN106097278A CN106097278A (zh) 2016-11-09
CN106097278B true CN106097278B (zh) 2021-11-30

Family

ID=57252842

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610474921.0A Active CN106097278B (zh) 2016-06-24 2016-06-24 一种多维信号的稀疏模型、重建方法和字典训练方法

Country Status (1)

Country Link
CN (1) CN106097278B (zh)

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107067380B (zh) * 2017-03-28 2020-04-28 天津大学 基于低秩张量和层次化字典学习的高分辨率图像重构方法
CN107480720B (zh) * 2017-08-18 2021-02-26 成都通甲优博科技有限责任公司 人体姿态模型训练方法及装置
CN107561576B (zh) * 2017-08-31 2023-10-20 中油奥博(成都)科技有限公司 基于字典学习正则化稀疏表示的地震信号恢复方法
CN108305297A (zh) * 2017-12-22 2018-07-20 上海交通大学 一种基于多维张量字典学习算法的图像处理方法
CN111212198B (zh) * 2018-11-22 2022-04-08 中移（杭州）信息技术有限公司 一种视频去噪的方法和设备
CN109921799B (zh) * 2019-02-20 2023-03-31 重庆邮电大学 一种基于聚能量字典学习的张量压缩方法
CN109934178A (zh) * 2019-03-18 2019-06-25 电子科技大学 一种基于Kronecker基稀疏表示的红外弱小目标检测方法
CN110969638B (zh) * 2019-11-12 2023-09-29 桂林电子科技大学 一种基于张量的背景减除方法及系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103279959A (zh) * 2013-06-13 2013-09-04 北京工业大学 一种二维分析稀疏模型、其字典训练方法和图像去噪方法
CN103279932A (zh) * 2013-06-05 2013-09-04 北京工业大学 一种二维合成稀疏模型和基于该模型的字典训练方法
CN105046672A (zh) * 2015-06-30 2015-11-11 北京工业大学 一种图像超分辨率重建方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103279932A (zh) * 2013-06-05 2013-09-04 北京工业大学 一种二维合成稀疏模型和基于该模型的字典训练方法
CN103279959A (zh) * 2013-06-13 2013-09-04 北京工业大学 一种二维分析稀疏模型、其字典训练方法和图像去噪方法
CN105046672A (zh) * 2015-06-30 2015-11-11 北京工业大学 一种图像超分辨率重建方法

Also Published As

Publication number Publication date
CN106097278A (zh) 2016-11-09

Similar Documents

Publication Publication Date Title
CN106097278B (zh) 一种多维信号的稀疏模型、重建方法和字典训练方法
Yokota et al. Smooth PARAFAC decomposition for tensor completion
US11221990B2 (en) Ultra-high compression of images based on deep learning
Patel et al. Kernel sparse subspace clustering
Figueiredo et al. Majorization–minimization algorithms for wavelet-based image restoration
Qi et al. Multi-dimensional sparse models
CN107067367A (zh) 一种图像超分辨重建处理方法
Qi et al. Using the kernel trick in compressive sensing: Accurate signal recovery from fewer measurements
CN112699899A (zh) 一种基于生成对抗网络的高光谱图像特征提取方法
CN111640059B (zh) 基于高斯混合模型的多字典图像超分辨方法
CN113205449A (zh) 表情迁移模型的训练方法及装置、表情迁移方法及装置
Yger et al. Supervised logeuclidean metric learning for symmetric positive definite matrices
CN105513033A (zh) 一种非局部联合稀疏表示的超分辨率重建方法
Qi et al. Two dimensional synthesis sparse model
CN110097499B (zh) 基于谱混合核高斯过程回归的单帧图像超分辨率重建方法
Wen et al. The power of complementary regularizers: Image recovery via transform learning and low-rank modeling
Silva et al. Efficient separable filter estimation using rank-1 convolutional dictionary learning
Zhang et al. Tensor Recovery Based on A Novel Non-convex Function Minimax Logarithmic Concave Penalty Function
Yao A compressed deep convolutional neural networks for face recognition
Liu et al. Hyperspectral image classification based on extended morphological profile features and ghost module
CN113554569A (zh) 基于双记忆字典的人脸图像复原系统
Chen et al. Tensor network for image classification
CN112966753A (zh) 一种基于熵稳定约束的自动数据降维方法
Luo et al. Convolutional Laplacian sparse coding
CN111475768A (zh) 一种基于低相干单位范数紧框架的观测矩阵构造方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant