基于信息熵稀疏的深度非负矩阵解混方法和存储介质

1.本发明涉及图像处理技术技术领域，具体涉及一种基于信息熵稀疏的深度非负矩阵解混方法和存储介质。


背景技术：

2.高光谱解混技术是一项重要的图像处理技术，被广泛应用于农业、矿产勘探和环境监测等各领域。其中，非负矩阵分解(nmf)算法是研究的热门，但nmf算法有明显的非凸性，存在多个局部极小值。且大多数基于nmf的解混方法只考虑单层中的信息，而忽略了原始数据包含具有隐藏信息的分层特征的事实，对于复杂和高度混合的数据而言，这可能会导致解混性能不理想。


技术实现要素：

3.本发明解决的一个主要问题是现有技术中对于复杂和高度混合的数据而言，基于nmf的解混方法解混性能不理想的问题。
4.根据本发明的一个方面，本发明提供一种基于信息熵稀疏的深度非负矩阵解混方法，包括：
5.对端元矩阵和丰度矩阵进行vca
‑
fcls初始化，并建立深度nmf解混模型；
6.将初始化后的所述端元矩阵和初始化后的所述丰度矩阵作为第一目标函数的输入，并利用所述第一目标函数，基于乘性迭代规则对所述深度nmf解混模型预训练，所述第一目标函数为：
[0007][0008]
其中，a为光谱特征矩阵,s为丰度矩阵，x为具有b个波段和p个像元的高光谱数据集，s
ij
代表第i种地物类型在第j个像元的丰度，i＝1,2,...,m，j＝1,2,...,p；λ是正则化参数；
[0009]
利用梯度下降法迭代更新预训练后的所述端元矩阵和所述丰度矩阵；
[0010]
将更新后的所述端元矩阵和更新后的所述丰度矩阵作为第二目标函数的输入，并利用所述第二目标函数对预训练后的所述nmf解混模型迭代微调，所述第二目标函数为：
[0011][0012]
其中，minf(a,s)表示重构误差，l表示第l层，l表示最后一层，代表第l层的重构的丰度矩阵，ψ
l
‑1＝a1a2…
a
l
‑1，且对于ψ
l
‑1，当且仅当l＝1时，代表一个单位矩阵；
[0013]
再次迭代更新微调后的所述端元矩阵和所述丰度矩阵，并通过数据分析确定所述第二目标函数的最优参数；
[0014]
根据所述最优参数获得所述nmf解混模型的全局最优解。
[0015]
进一步地，所述深度nmf解混模型的数学表达式为：
[0016]
x≈a1a2…
a
l
‑1a
l
s
l
。
[0017]
进一步地，所述利用所述第二目标函数对预训练后的所述nmf解混模型迭代微调还包括：
[0018]
对所述深度nmf结构的每个中间层结果均进行微调；
[0019]
进一步地，所述基于乘性迭代规则对所述深度nmf解混模型预训练还包括：
[0020]
重新定义高光谱数据集为第一表达式：
[0021][0022]
其中，δ为一常量，为新的高光谱数据集，向量1是一个全为1的行向量；定义光谱特征矩阵为第二表达式：
[0023][0024]
为新的光谱特征矩阵。
[0025]
进一步地，所述基于乘性迭代规则对所述深度nmf解混模型预训练还包括：
[0026]
基于所述第一目标函数分别对所述光谱特征矩阵和所述丰度矩阵求偏导，
[0027]
并利用梯度下降法对所述光谱特征矩阵和所述丰度矩阵进行迭代，得到所述光谱特征矩阵和所述丰度矩阵的乘性迭代规则。
[0028]
进一步地，所述基于乘性迭代规则对所述深度nmf解混模型预训练还包括：
[0029]
基于所述第二目标函数分别对第l层所述光谱特征矩阵和第l层所述丰度矩阵求偏导，利用梯度下降法进行迭代，得到所述深度nmf解混模型微调第l层分解结果的乘性迭代公式。
[0030]
进一步地，所述利用所述第二目标函数对预训练后的所述nmf解混模型迭代微调还包括：
[0031]
微调的迭代过程利用整个所述深度nmf解混模型的参数，并且以预训练得到的所述端元矩阵和所述丰度矩阵每层的结果为迭代的初始值。
[0032]
根据本发明的另一个方面，还提供一种存储介质，所述存储介质为计算机可读存储介质，所述计算机可读存储介质上存储有如前任一所述的一种基于信息熵稀疏的深度非负矩阵解混方法。
[0033]
本发明通过使用深度nmf模型代替传统的单层nmf模型，并进行信息熵稀疏，提出了一种基于信息熵稀疏的深度nmf算法(hdnmf)，能够更好的挖掘隐含层信息，充分考虑空间信息，通过模拟数据和真实cuprite数据集进行验证，验证了该算法的有效性。
附图说明
[0034]
本发明构成说明书的一部分附图描述了本发明的实施例，并且连同说明书一起用于解释本发明的原理。
[0035]
图1为本发明实施例中深度nmf解混模型示意图。
[0036]
图2为本发明实施例中dnmf算法流程图示意图。
[0037]
图3为本发明实施例中模拟数据4种端元的光谱曲线示意图。
[0038]
图4为本发明实施例中合成高光谱数据示意图。
[0039]
图5为本发明实施例中模拟数据端元丰度对比图。
[0040]
图6为本发明实施例中变化时sad和rmse值变化图。
[0041]
图7为本发明实施例中hdnmf在不同层数下的sad和rmse。
[0042]
图8为本发明实施例中真实数据中玉髓和榍石的丰度图像。
具体实施方式
[0043]
下面将结合附图来详细描述本发明的各种示例性实施例。应注意到：除非另外具体说明，否则在这些实施例中阐述的部件和步骤的相对布置、数字表达式和数值不限制本发明的范围。
[0044]
同时，应当明白，为了便于描述，附图中所示出的各个部分的尺寸并不是按照实际的比例关系绘制的。
[0045]
以下对至少一个示例性实施例的描述实际上仅仅是说明性的，决不作为对本发明及其应用或使用的任何限制。
[0046]
为使本发明的目的、技术方案和优点更加清楚明白，以下结合具体实施例，并参照附图，对本发明进一步详细说明。
[0047]
对于相关领域普通技术人员已知的技术、方法和设备可能不作详细讨论，但在适当情况下，所述技术、方法和设备应当被视为授权说明书的一部分。
[0048]
在这里示出和讨论的所有示例中，任何具体值应被解释为仅仅是示例性的，而不是作为限制。因此，示例性实施例的其它示例可以具有不同的值。
[0049]
应注意到：相似的标号和字母在下面的附图中表示类似项，因此，一旦某一项在一个附图中被定义，则在随后的附图中不需要对其进行进一步讨论。
[0050]
实施例一，如图1所示，为本实施例中的深度nmf模型结构示意图，常用的光谱混合模型有线性混合模型和非线性混合模型两种，本发明实施例一使用的是线性混合模型，它具有清晰的物理特性，有较好的理论基础.假设像元的光谱曲线是多个端元光谱特征的线性组合，按照各个端元的丰度系数加权得到的.设为具有b个波段和p个像元的高光谱数据集，则该数据集可以用线性混合模型表示为：
[0051]
x＝as+e
ꢀꢀ
(1)
[0052]
其中，代表所有地物类型的光谱特征矩阵，具有b个波段和m个端元，并且a＞0；s
ij
代表第i种地物类型在第j个像元的丰度，为了让解混的结果具有物理意义，s中每一列的m个元素之和都为1；e代表传感器接收到的噪声干扰。
[0053]
lee和seung根据噪声服从不同的概论分布提出了基于欧式距离和基于k
‑
l散度的求解方法，
[0054]
min f(a,s)＝||x
‑
as||
f s.t.a≥0,s≥0
ꢀꢀ
(2)，
[0055]
其中，端元矩阵a、丰度矩阵s的所有元素应同时满足非负约束(anc)，并且在高光谱图像解混中，丰度矩阵s应额外满足列向量和为一约束(asc)。
[0056]
由于nmf算法具有明显的非凸性，存在多个局部极小值，需要进一步优化。目前，大多的nmf算法都是通过稀疏正则化来进行优化。在本算法中，为了使最终结果可以接近全局最优解，通过顶点成分分析(vca)来初始化端元矩阵，用完全约束最小二乘法(fcls)来初始化丰度矩阵，提高迭代的收敛速度。
[0057]
利用信息熵来对丰度矩阵进行稀疏约束，同时借助深度学习的思想，将单层的nmf模型扩展到深度nmf模型。如图1所示，是深度nmf模型的框架图。从中可以看出，它的结构是由多个单层的nmf模型一层接着一层堆叠而成的。对于第l层，将前一层经过非负矩阵分解后的输出s
l
‑1作为当前nmf层的输入，分解得到a
l
和s
l
，以此类推，直到最后一层第l层。
[0058]
那么，可以推导出深度nmf模型的数学表达式为：
[0059]
x≈a1a2…
a
l
‑1a
l
s
l
ꢀꢀ
(3),
[0060]
本实施例当中的深度nmf模型解混方法主要由两个阶段组成，分别是预训练阶段和微调阶段。
[0061]
在预训练阶段，其目标函数(即第一目标函数)如下：
[0062][0063]
其中，min f(a,s)表示重构误差，a为光谱特征矩阵，s为丰度矩阵，x为具有b个波段和p个像元的高光谱数据集，s
ij
代表第i种地物类型在第j个像元的丰度，i＝1,2,...,m，j＝1,2,...,p；λ是正则化参数。
[0064]
由于预训练只能保证每个单层的nmf是收敛的，但是最后得到的解混结果是对整体进行的计算，这样可能会导致误差的累积，为了降低深度nmf结构的整体误差，还需要通过某种优化法则分别对每个中间层结果进行微调。
[0065]
在微调阶段，将公式(3)中的||x
‑
as||重写为更详细的形式，如下：
[0066][0067]
为了方便理解和后续的求解，需要定义两个新的变量：
[0068]
ψ
l
‑1＝a1a2…
a
l
‑1ꢀꢀ
(6)
[0069][0070]
其中，代表第l层的重构的丰度矩阵，ψ
l
‑1＝a1a2…
a
l
‑1，且对于ψ
l
‑1，当且仅当l＝1时，代表一个单位矩阵；
[0071]
因此，在微调阶段，目标函数如下(即第二目标函数)：
[0072][0073]
其中，l表示第l层，l表示最后一层，代表第l层的重构的丰度矩阵，ψ
l
‑1＝a1a2…
a
l
‑1，且对于ψ
l
‑1，当且仅当l＝1时，代表一个单位矩阵。
[0074]
在高光谱图像解混中，为了让丰度矩阵s满足asc约束，需要重新定义数据集和光谱特征矩阵为：
[0075][0076]
其中δ控制和为一约束的影响，向量1是一个全为1的行向量，因此，可以得到a和s的乘性迭代规则。
[0077]
在预训练阶段，第一目标函数分别对a，s求偏导得：
[0078][0079][0080]
其中，s
t
是矩阵s的转置，a
t
是矩阵a的转置，*是指矩阵的元素对应相乘；选择合适的步长，利用梯度下降法进行迭代,得到a和s的乘性迭代规则如下:
[0081]
a＝a*xs
t
/(ass
t
+ε)
ꢀꢀ
(12)
[0082][0083]
式中利用小正数ε来使分数恒为正数，当迭代到一定次数，f(a，s)的值趋于稳定，其中为新的光谱特征矩阵，是的转置。
[0084]
在微调阶段，对微调阶段的目标函数(8)分别对a
l
、s
l
求偏导，选择合适的步长，例如:a
←
a
‑
步长*对a的偏导，那么会使步长为a/asst，这样相乘时抵消了a，a
←
a*xst/asst；
[0085]
利用梯度下降法对a
l
、s
l
进行迭代，推导得到深度nmf微调第l层分解结果的乘性迭代公式为：
[0086][0087][0088]
其中，是的转置，a
t
是矩阵a的转置，*是指矩阵的元素对应相乘。
[0089]
从式子(15)可以看出，和预训练不一样的是：微调的迭代过程是以预训练得到的结果为初始值，并且用到了整个深度nmf模型的参数，当前层参数的更新不仅跟它之前的所有层相关，还与其之后的所有层相关，这就保证了每一次迭代更新后的分解结果都是基于整体误差推算出来的，降低了误差的累积。
[0090]
该算法的整体流程如图2所示：
[0091]
其中，在预训练和微调阶段的迭代次数选择为1000次，正则参数λ从有限集合{0.0001，0.0005，0.001，0.005，0.01，0.05，0.1，0.2}中选择。
[0092]
实施例一的实验验证：
[0093]
1.模拟数据
[0094]
模拟数据由4种从美国地质调查局(usgs)数字光谱库中挑选的4种端元光谱线性混合而成，光谱值表示相应地物在0.4
‑
2.5μm波长范围内的反射率，波段数均为b＝224，本文只选择了188个高信噪比的波段用于实验，图3展示了4种端元的光谱特征。
[0095]
合成影像包含48
×
48个像元，每个像元包含188个波段。合成影像的丰度图像满足anc以及asc，其中第1行为纯像元，第2至4行分别为2种、3种、4种端元的混合，方块周围的背
景像素也为相同的4种端元以不同的比例的混合，模拟数据的合成影像如图4所示：
[0096]
2.真实数据
[0097]
在本文中，使用的真实数据是1997年6月在美国内华达州拍摄的cuprite高光谱数据，平均有224个波段。在本文中截取了大小为143
×
123像素的遥感数据，选择了188个高信噪比的波段用于实验。对该图像进行mnf变换和pca变换，确定该区域主要包含12种地物，即端元数m＝12。
[0098]
3.实验性能指标
[0099]
对于模拟数据的解混效果的评估，由于我们知道模拟数据所包含端元的光谱特征和丰度系数等先验信息，所以可以分别使用光谱角距离(sad)和均方根误差(rmse)来定量评估提取的端元和丰度的精度，具体定义如下：
[0100]
(1)光谱角距离sad：度量原始端元光谱特征a与其估计值之间的差异，公式如下：
[0101][0102]
光谱角距离用来评价实验所得的端元光谱矩阵与参考值的近似度。由上述定义可知：sad值越小，与参考结果越相似，实验结果越好。
[0103]
(2)均方根误差rmse：度量丰度真值s
ij
与丰度估计值之间的差异，公式如下：
[0104][0105]
均方根误差用来评估实验所得的丰度矩阵与参考值的近似度。由上述定义可知：rmse值越小，与参考结果越相似，实验结果越好。
[0106]
对于真实数据的解混效果的评估，一般通过定性的方法进行评估。目视评估解混结果的好坏，观察提取出来端元的丰度的物理意义是否明确，丰度构成是否存在明显偏差，以及观察解混得到的光谱特征曲线是否符合某种地物的特征。
[0107]
实验结果与分析
[0108]
1.模拟数据实验
[0109]
本实验首先用合成数据集来验证hdnmf算法的有效性，并与chnmf算法和l1/2
‑
nmf算法进行比较。实验中的所有结果取重复20次的均值。
[0110]
如图4所示，是合成数据的4种端元通过l1/2
‑
nmf，chnmf和hdnmf这三种算法解混出来的丰度图与参考丰度图的比较。端元丰度图中颜色的深浅分别表示该区域含相应端元的比例，即丰度值的大小，颜色越深表示所含端元比例越少(丰度值小)，深色表示不含该端元，亮色(浅色)表示此区域全部都是该种端元，不含其他端元。由图5可得，hdnmf算法的解混效果最好。
[0111]
表1显示了这三种算法解混结果的sad值和emse值的比较。从实验结果上看，hdnmf算法的两个评估指标均要优于l1/2
‑
nmf算法和chnmf算法。
[0112]
表1三种nmf解混算法在模拟数据中的sad值和rmse值
[0113][0114]
2.参数的敏感性分析
[0115]
本发明提出的hdnmf算法涉及到参数λ的变化，首先设定深度层数为3，参数λ从有限集合{0.0001，0.0005，0.001，0.005，0.01，0.05，0.1，0.2}中选择，随着λ的变化，sad和rmse的变化情况分别如图6中的图6a和图6b所示。由图6观察可得，当λ＝0.001时，各地物都有较好的解混效果。
[0116]
3.深度层数的影响分析
[0117]
深度学习层数是hdnmf算法的重要影响因素，本实验中，设定λ＝0.001，层数设置从1到5，各个层的sad和rmse值如图7所示。当层数为3时，rmse和sad的值都最小，解混效果最好，并且随着层数的增加，解混精度反而变差，这是因为随着学习层数的增加，累积了误差，导致了解混效果的降低。因此，在本文的实验中选取3层结构来进行实验。
[0118]
4.真实数据实验
[0119]
通过真实数据实验，进一步对比分析l1/2
‑
nmf，chnmf和hdnmf这三种算法的解混情况。从实验结果中选取了玉髓和榍石这两种比较典型的地物类型用来做解混算法的对比。
[0120]
表2三种解混算法在真实数据中的sad值
[0121][0122]
如图8所示，给出了l1/2
‑
nmf、chnmf、hdnmf三种算法提取的两种端元丰度图。对端元估计值和usgs参考光谱进行比较，结果如表3所示。可以看出，hdnmf算法与l1/2
‑
nmf、chnmf算法相比，hdnmf算法的解混效果最好。
[0123]
以上所述仅为本发明的较佳实施例而已，并不用于限制本发明，凡在本发明的精神和原则范围之内所作的任何修改、等同替换以及改进等，均应包含在本发明的保护范围
之内。
[0124]
还需要说明的是，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个
……”
限定的要素，并不排除在包括要素的过程、方法、商品或者设备中还存在另外的相同要素。