遥感技术自20世纪以来取得了巨大进展，高光谱遥感技术作为其中之一已经广泛应用到各个领域。高光谱成像技术将成像技术与光谱技术相结合，可以同时获得空间连续和光谱连续的数据，是一种具有高维特征空间的非线性结构数据，并且提供了丰富的光谱信息。传统高光谱分类任务利用光谱信息进行分类，使用的分类算法包括平行六面体分类、K近邻、最大似然、最小距离和逻辑回归(Foody和Mathur，2004)。目前虽然有许多针对高光谱图像分类的研究(冉琼等，2018；方帅等，2019)，但是由于高光谱图像标签数据少，标注困难，使得无监督聚类任务成为高光谱图像研究的重要课题之一(于文博等，2016；Landgrebe，2002)。

虽然已经提出了大量的聚类方法(Kohonen，1990；Reynolds，2015；Ester等，1996；Zhang等，2017a；Zhang等，2017b), 如K均值聚类算法(K-means)(Hartigan和Wong，1979)和谱聚类算法(SC)(von Luxburg，2007)，但是由于这些传统方法作为浅层线性模型，使用的相似性测度效率较低，对于具有高维特性或者具有非线性结构特性的高光谱数据，其聚类效果较差。传统聚类算法大多数存在“维数灾难”问题。研究者们针对上述问题提出了高维数据聚类方法，如子空间聚类；多种特征提取方法，如主成分分析(PCA)、独立成分分析(ICA)和小波变换等。

子空间聚类方法的目的是为未标记数据集寻找一组低维子空间，从而将其划分为不同的簇(Vidal, 2011)。如稀疏子空间聚类算法(SSC)(Elhamifar和Vidal，2009)和低秩子空间聚类算法(LRSC)(Liu等，2012)，SSC算法是SC算法框架下的稀疏优化求解，近年来许多基于稀疏表示和低秩近似的子空间聚类方法得到了广泛关注。这些方法的关键是找到数据的稀疏和低秩表示，然后在稀疏系数矩阵上建立相似图。虽然一些子空间聚类方法(Nie等，2011；Yang等，2016；Li和Vidal，2016；Peng等，2017；Yuan等，2016；Li等，2018)取得了显著效果，但是存在如下局限：首先这些子空间聚类方法的聚类效果依赖于亲和矩阵；其次由于这些方法对数据做了线性假设，导致它们无法有效处理非线性结构的数据。特征提取方法PCA(Wold等，1987)虽然可以有效地降维并且提取光谱特征，但是作为线性变换方法，所提取的光谱特征的有效性仍然面临挑战。为了提取更加有效的特征，研究者们针对高光谱数据非线性结构的特点提出了一些核方法：如核子空间聚类方法，核PCA(Hofmann等，2008)等。它们首先将数据映射到预先定义的核空间中，然后在核空间中执行子空间聚类或特征提取。这些核方法的一个主要缺点是它们的性能严重依赖于使用的核函数。然而，面对不同类型的核函数，选择哪种核函数主要依赖于经验，并且核函数所定义的隐式特征空间是否适用于聚类并没有明确的理论依据。

随着深度学习(Schmidhuber，2015)的发展, 深度神经网络(DNN)可以将数据非线性映射到新的特征空间中，在提取深层特征和降维方面具有传统聚类算法无可比拟的优势。深度自编码器网络(DAE)作为典型的无监督学习人工神经网络，具有对称的网络结构，是无监督表示学习中最重要的算法之一。它是一种功能强大的映射函数训练方法，保证了重构数据和原始数据之间的重构误差最小。由于隐藏层的维数通常小于数据层，因此它可以帮助发现数据最显著的特征。虽然在有监督学习中，DAE主要用于寻找更好的参数初始化，但将其与无监督聚类相结合也是很自然的。基于DAE的聚类算法(Chang等，2017)是一种流行的深度聚类算法模型，如基于自编码器的数据聚类(AEKM)(Song等，2013)。深度聚类算法模型通常分为两步：1)通过自编码器的编码过程对数据进行降维和特征提取，得到低维特征空间；2)在低维特征空间中对数据的低维本质特征进行聚类。然而以往的工作将特征提取和聚类任务剥离开，二者独立进行。由于特征提取过程缺乏与聚类任务相关的约束，这些方法仅仅基于重构误差学习数据深层表示，得到的数据特征适合重构原始图片，而未必适合聚类任务。并且深度聚类算法多应用于图像聚类，如人脸数据集、自然图像数据集(手写数字数据集、物体图像数据集等)，但是在高光谱图像数据集上的应用很少。

为了解决上述问题，充分挖掘高光谱图像的光谱信息，本文提出了一种嵌入式深度神经网络模糊C均值聚类(EDFCC)算法。所提方法能够有效挖掘高光谱图像的光谱信息，并用于高光谱图像聚类。EDFCC算法将模糊C均值聚类算法(FCM)嵌入到DAE网络结构中，同时优化特征提取过程和聚类过程。首先，在特征提取过程中，为了保持局部结构，通过最小化重构误差学习高光谱图像每个像素点的特征表示，即通过DAE网络完成光谱特征提取过程；其次，为了自适应学习适用于聚类的光谱特征，更好地对高光谱图像进行聚类，以FCM聚类算法约束特征提取过程。EDFCC算法具有两个任务联合优化的优势，并且利用DAE网络强大的降维和非线性映射能力将原始数据映射到潜在特征空间，同时以聚类误差约束特征映射过程，从而得到更加适合聚类的深层特征。

本文提出的EDFCC算法有如下贡献：1)EDFCC算法通过DAE网络学习原始数据的深层表示，其特征映射过程是非线性的，能够处理非线性结构的高维数据。2)利用聚类误差约束DAE网络的编码过程，并联合重构误差对网络参数进行优化，从而使学到的特征更加适合于聚类。3)动态调整聚类指示矩阵，并且不需要额外的训练过程，提高了聚类效率。

一些子空间聚类任务需要将数据从高维特征空间降至低维特征空间，但是当前的子空间聚类算法大多是线性模型(Patel和Vidal，2014；Yin等，2016；Elhamifar和Vidal，2013；Peng等，2015)，在很多情况下，无法对非线性结构数据进行建模。由于深度神经网络具有强大的非线性建模和特征表示能力，近年来提出了一些基于深度神经网络的聚类方法。然而，当前深度聚类算法的特征提取和聚类是两个相对独立的过程，聚类算法并未嵌入到特征提取过程中，因此，聚类算法无法约束特征提取过程。针对这一问题，本文提出了EDFCC算法，该方法以DAE网络为基础，将FCM聚类算法嵌入至神经网络中，约束EDFCC算法的特征提取过程。

深度学习以强大的学习表征能力在众多应用中取得了巨大成功。在无监督学习领域，良好的数据特征和数据表达尤为重要。基于此，提出了无监督深度神经网络和聚类联合学习的模型，希望更好地解决无监督聚类问题。本小节主要阐述使用DAE网络提取高光谱图像深层光谱特征的原因以及基础DAE(Hinton和Salakhutdinov，2006)网络。

一幅典型的高光谱图像的场景覆盖范围广、通道数量多、往往在2~3位数量级。高光谱图像的聚类任务是根据像素点各自的上百个光谱通道数据为每个像素点赋予一个类别标签，所以高光谱图像的聚类任务是基于像素的。

引入DAE网络提取高光谱图像深层光谱特征，的主要原因是：1)由于光照条件复杂多变，同一类地物在不同位置表现出不同的光谱特征。例如，同一块草坪直接暴露在太阳下的光谱特征与被高层建筑遮挡的光谱特征不同。2)不同地物在各自的光照条件下可能表现出相同或相似的光谱特征。3)除了自身的散射，地物的光谱特征还会受到周围环境散射的影响。4)地物的光谱特征还会受到其他因素的影响，如传感器误差、大气散射等。

考虑到上述因素，一个特定类别的地物在特征空间中存在多个方向的变化，在复杂的实际情况下，这些复杂的光谱变化使高鲁棒性和不变性的深层特征尤为重要。

深度神经网络提取到的深层特征是其浅层特征的组合和抽象，而抽象的特征大多数对输入具有局部不变性。因此，为了获得具有不变性的深层特征，并解决上述问题，针对聚类任务使用DAE网络提取高光谱图像光谱特征的深层表示。

DAE网络编码过程和解码过程都不涉及标签信息，因此它是一种无监督的神经网络模型，由于它具有良好的特征学习和特征提取能力，常用于无监督学习领域作为特征提取器。自编码器具有对称的网络结构，一个自编码器包括编码层和解码层，并且编码层的层数与解码层的层数是相等的。编码层的作用是将原始数据非线性映射到一个隐藏空间，由于隐藏层的维数通常小于输入层，因此它可以帮助找到数据最显著的特征，得到数据的低维表示，解码层的作用是根据学习到的数据低维表示重构出相应的原始输入数据。

DAE网络共由$M$+1层叠加而成，其中包括1个输入层和由$M$个非线性变换得到的$M$个输出层，DAE网络以高光谱数据像元光谱向量作为输入样本，即输入层定义为$\boldsymbol{h}^{(0)} \in {\bf{R}}^{d}$，共执行$M$次非线性映射，$\boldsymbol{h}^{(0)}$即初始像元光谱向量，那么DAE后续层的定义为

$ {\mathit{\boldsymbol{h}}^{\left( m \right)}} = g\left( {{\mathit{\boldsymbol{W}}^{\left( m \right)}}{\mathit{\boldsymbol{h}}^{\left( {m - 1} \right)}} + {\mathit{\boldsymbol{b}}^{\left( m \right)}}} \right) \in {{\bf{R}}^{{d^{\left( m \right)}}}} $

(1)

式中，$\boldsymbol{h}^{(m)}$表示DAE网络第$m$层的输出($m$=1, 2, …, $M$为层数)，$g$(·)表示非线性激活函数，${d^{(m)}}$表示第$m$层输出数据所具有的维度，${\mathit{\boldsymbol{W}}^{(m)}} \in {{\bf{R}}^{{d^{(m)}} \times {d^{(m - 1)}}}}$表示第$m$层的权重矩阵，${\mathit{\boldsymbol{b}}^{(m)}} \in {{\bf{R}}^{{d^{(m)}}}}$表示第$m$层的偏置。特别地，如果给定$\boldsymbol{h}^{(0)}$作为网络的输入层，那么网络最后一层的输出为

$ {\mathit{\boldsymbol{h}}^{\left( M \right)}} = g\left( {{\mathit{\boldsymbol{W}}^{\left( M \right)}}{\mathit{\boldsymbol{h}}^{\left( {M - 1} \right)}} + {\mathit{\boldsymbol{b}}^{\left( M \right)}}} \right) \in {{\bf{R}}^{{d^{\left( M \right)}}}} $

(2)

通过式(1)可以发现，整个DAE网络是一种由权重、偏置和激活函数决定的非线性映射, 并且${{\bf{R}}^d} = {{\bf{R}}^{{d^{(M)}}}}$。如果有$N$个样本点，$N$表示高光谱图像的像元数量，则定义输入的样本集合、中间层编码集合(深层光谱特征集合)和输出层重构集合分别为

$ {\mathit{\boldsymbol{H}}^{\left( 0 \right)}} = \left[ {\mathit{\boldsymbol{h}}_1^{\left( 0 \right)},\mathit{\boldsymbol{h}}_2^{\left( 0 \right)}, \cdots ,\mathit{\boldsymbol{h}}_N^{\left( 0 \right)}} \right] \in {{\bf{R}}^{N \times d}} $

(3)

$ {\mathit{\boldsymbol{H}}^{\left( {\frac{M}{2}} \right)}} = \left[ {\mathit{\boldsymbol{h}}_1^{\left( {\frac{M}{2}} \right)},\mathit{\boldsymbol{h}}_2^{\left( {\frac{M}{2}} \right)}, \cdots ,\mathit{\boldsymbol{h}}_N^{\left( {\frac{M}{2}} \right)}} \right] \in {{\bf{R}}^{N \times d\left( {\frac{M}{2}} \right)}} $

(4)

$ {\mathit{\boldsymbol{H}}^{\left( M \right)}} = \left[ {\mathit{\boldsymbol{h}}_1^{\left( M \right)},\mathit{\boldsymbol{h}}_2^{\left( M \right)}, \cdots ,\mathit{\boldsymbol{h}}_N^{\left( M \right)}} \right] \in {{\bf{R}}^{N \times d\left( M \right)}} $

(5)

DAE网络结构如图 1所示，其中${\mathit{\boldsymbol{H}}^{(0)}}$表示原始高光谱图像，${\mathit{\boldsymbol{H}}^{(M)}}$表示重构的高光谱图像，${\mathit{\boldsymbol{H}}^{\left({\frac{M}{2}} \right)}}$表示网络模型学习到的高光谱图像深层光谱特征。

基础的DAE网络有3种结构，分别是压缩结构、稀疏结构和等效维结构。当隐藏层神经元数目小于输入层神经元数目时，称其为压缩结构；当隐藏层神经元数目大于输入层神经元数目时，称其为稀疏结构；当隐藏层神经元数目等于输入层神经元数目时，称其为等效维结构。其损失函数, 即图 1中的重构误差可以是均方误差或者是交叉熵误差，分别定义为

$ \begin{array}{*{20}{c}} {{J_{{\rm{rec}}}}\left( {\mathit{\boldsymbol{W}},\mathit{\boldsymbol{b}}} \right) = }\\ {\sum\limits_{i = 1}^n J \left( {\mathit{\boldsymbol{h}}_i^{\left( 0 \right)},\mathit{\boldsymbol{h}}_i^{\left( {M} \right)}} \right) = \sum\limits_{i = 1}^n {\left\| {\mathit{\boldsymbol{h}}_i^{\left( 0 \right)} - \mathit{\boldsymbol{h}}_i^{\left( {M} \right)}} \right\|_2^2} } \end{array} $

(6)

$ \begin{array}{*{20}{c}} {{J_{\rm{c}}}\left( {\mathit{\boldsymbol{W}},\mathit{\boldsymbol{b}}} \right) = \sum\limits_{i = 1}^n J \left( {\mathit{\boldsymbol{h}}_i^{\left( 0 \right)},\mathit{\boldsymbol{h}}_i^{\left( {M} \right)}} \right) = }\\ { - \sum\limits_{i = 1}^n {\left( {\mathit{\boldsymbol{h}}_i^{\left( 0 \right)}{{\log }_2}\left( {\mathit{\boldsymbol{h}}_i^{\left( {M} \right)}} \right) + \left( {1 - \mathit{\boldsymbol{h}}_i^{\left( 0 \right)}} \right){{\log }_2}\left( {1 - \mathit{\boldsymbol{h}}_i^{\left( {M} \right)}} \right)} \right)} } \end{array} $

(7)

k-means聚类算法在实现过程中，将每个样本划分到单一的类别中，即每个样本只能属于一种类别，这样的划分称为硬划分。与k-means聚类算法这种硬划分聚类算法不同，FCM是一种软划分聚类算法。它通过优化模糊目标函数得到每个样本对聚类中心的隶属度，从而决定样本点的归属。

FCM聚类算法允许单个样本点属于两个甚至更多的簇。该聚类算法由Dunn(1973)提出，并多用于模式识别领域(Bezdek，1981)。在本文提出的EDFCC算法中，FCM聚类算法以DAE网络模型学习到的高光谱图像深层光谱特征${\mathit{\boldsymbol{H}}^{\left({\frac{M}{2}} \right)}}$作为输入，以将函数

$ \begin{array}{l} {J_{\rm{f}}} = \sum\limits_{i = 1}^N {\sum\limits_{j = 1}^K {u_{ij}^m\left\| {\mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)} - {\mathit{\boldsymbol{c}}_j}} \right\|_2^2} } ,\;\;\;\;\;1 < m < \infty \\ {\rm{s}}.\;{\rm{t}}.\;\;\;\;\sum\limits_{j = 1}^K {{u_{ij}}} = 1 \end{array} $

(8)

最小化为目标。式(8)中，$m$表示大于1的实数，$N$表示高光谱图像的像元数量，$K$表示类别数，$\boldsymbol{h}_{i}^{\left(\frac{M}{2}\right)} \in {\bf{R}}^{d\left(\frac{M}{2}\right)}$表示高光谱图像第$i$个像元的深层光谱特征，${\mathit{\boldsymbol{c}}_j} \in {{\bf{R}}^d}^{\left({\frac{M}{2}} \right)}$表示第$j$个聚类中心，${u_{ij}}$表示第$i$个像元的深层光谱特征${\mathit{\boldsymbol{h}}_i}^{\left({\frac{M}{2}} \right)}$对第$j$个簇的隶属度，$\left\| \cdot \right\|$表示测度深层光谱特征与聚类中心相似度的任意范数。聚类过程对上述目标函数进行迭代优化，对样本集进行模糊划分，通过

$ {u_{ij}} = \frac{1}{{\sum\limits_{k = 1}^c {{{\left( {\frac{{\left\| {\mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)} - {\mathit{\boldsymbol{c}}_j}} \right\|}}{{\left\| {\mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)} - {\mathit{\boldsymbol{c}}_k}} \right\|}}} \right)}^{\frac{2}{{m - 1}}}}} }} $

(9)

更新隶属度矩阵$\mathit{\boldsymbol{U}} \in {{\bf{R}}^{N \times K}}$, 并通过

$ {\mathit{\boldsymbol{c}}_j} = \frac{{\sum\limits_{i = 1}^N {u_{ij}^m} \cdot \mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)}}}{{\sum\limits_{i = 1}^N {u_{ij}^m} }} $

(10)

更新聚类中心矩阵$\mathit{\boldsymbol{C}} \in {{\bf{R}}^{K \times d}}^{\left({\frac{M}{2}} \right)}$。

EDFCC算法以DAE网络为基础，具体结构可以根据不同的场景而变化。

EDFCC算法的模型框架如图 2所示。实验中采用7层神经网络结构，第1层和最后1层的神经元数量为不同高光谱数据集的波段数，中间5层一律采用100、50、25、50、100个神经元的结构。

EDFCC算法实现的步骤为：首先利用DAE网络提取高光谱图像光谱信息的深层特征，在隐藏层上构造FCM聚类算法，使用深层特征作为聚类算法的输入，完成聚类任务。

具体地，将待聚类像素点的所有光谱值组成光谱向量作为EDFCC算法中DAE网络的输入数据，对每个像素点的光谱向量进行深层特征提取，并将其作为聚类算法和解码器的输入，同时得到重构误差和聚类误差，反向训练模型，直至模型收敛，得到最后的聚类结果。

该框架以DAE网络为基础，通过最小化重构误差(式(6))训练整个DAE网络，同时以聚类误差(式(8))约束DAE网络的特征提取过程，重构误差和聚类误差联合训练网络参数，从而得到更加适合聚类的光谱特征。根据上述定义，框架的目标函数可以定义为

$ \mathop {\min }\limits_{\mathit{\boldsymbol{W}},\mathit{\boldsymbol{b}}} J = {J_{{\rm{rec}}}} + {\lambda _1}{J_{\rm{f}}} + {\lambda _2}{J_{{\rm{reg}}}} $

(11)

式中，${J_{{\rm{rec}}}}$表示EDFCC算法中DAE网络模型得到的重构高光谱图像${\mathit{\boldsymbol{h}}_i}^{(M)}$与原始高光谱图像${\mathit{\boldsymbol{h}}_i}^{(0)}$之间的误差，其结果为标量；${J_{\rm{f}}}$表示以EDFCC算法中DAE网络模型中间层的深层光谱特征${\mathit{\boldsymbol{H}}^{\left({\frac{M}{2}} \right)}}$作为输入，在其上构造FCM聚类算法后预测类别与真实类别之间的误差，其结果为标量。${\lambda _2}$是正则化参数，${\lambda _1}$是用来平衡重构误差和聚类误差的参数，因为当重构误差很小的时候，得到的特征未必同样适用于聚类；同样当聚类误差很小的时候，得到的特征或许不能很好地重构出原始图片，所以${\lambda _1}$可以防止网络学习的参数过于偏重二者之一，${J_{{\rm{reg}}}}$为正则项，其结果为标量，表示为

$ {J_{{\rm{reg}}}} = \sum\limits_{i = 1}^M {\left( {\left\| {{\mathit{\boldsymbol{W}}_i}} \right\|_{\rm{F}}^2 + \left\| {{\mathit{\boldsymbol{b}}_i}} \right\|_2^2} \right)} $

(12)

网络正则是一种控制过拟合的手段，通过对网络模型参数添加正则项来约束网络参数，避免过拟合。具体地，通过使网络模型中每一层的权重矩阵$\mathit{\boldsymbol{W}}$和偏置$\mathit{\boldsymbol{b}}$的正则项(L₁范数，L₂范数，L₁ + L₂范数等)之和最小化，促使网络模型参数相对较小。由目标函数可知，若网络模型参数较大，数据稍有偏移就会大幅影响结果；反之若网络模型参数较小，即使数据变化范围比较大，对结果的干扰也较小，增强抗扰动能力，实现避免模型过拟合的目的，提高了模型的泛化能力。于是完整的目标函数定义为

$ \begin{array}{*{20}{c}} {\mathop {\min }\limits_{{\mathit{\boldsymbol{W}}_i},{\mathit{\boldsymbol{b}}_i}} J\left( {{\mathit{\boldsymbol{W}}_i},{\mathit{\boldsymbol{b}}_i}} \right) = \left\| {{\mathit{\boldsymbol{H}}^{(0)}} - {\mathit{\boldsymbol{H}}^{\left( M \right)}}} \right\|_{\rm{F}}^2 + }\\ {{\lambda _1}\left( {\sum\limits_{i = 1}^N {\sum\limits_{j = 1}^K {u_{ij}^m\left\| {\mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)} - {\mathit{\boldsymbol{c}}_j}} \right\|_2^2} } } \right) + }\\ {{\lambda _2}\sum\limits_{i = 1}^M {\left( {\left\| {{\mathit{\boldsymbol{W}}_i}} \right\|_{\rm{F}}^2 + \left\| {{\mathit{\boldsymbol{b}}_i}} \right\|_2^2} \right)} } \end{array} $

(13)

第1项${J_{{\rm{rec}}}}$通过最小化输入数据本身的重构误差保持局部性，输入数据本身充当了特征映射过程的监督者；第2项${J_{\rm{f}}}$通过最小化聚类误差约束非线性变换过程；第3项${J_{{\rm{reg}}}}$是正则项，避免网络模型过拟合。

更重要的是在EDFCC算法模型中，隶属度矩阵$\mathit{\boldsymbol{U}}$不是固定的，它会在模型的每一次迭代过程中重新计算，并且不需要额外的预处理过程，一步完成聚类，使得效率得以有效提升。

为了方便算法优化过程的推导，将式(13)重写为单样本形式，具体为

$ \begin{array}{*{20}{c}} {J = \frac{1}{2}\sum\limits_{i = 1}^n {\left( {\left\| {\mathit{\boldsymbol{h}}_i^{\left( 0 \right)} - \mathit{\boldsymbol{h}}_i^{\left( M \right)}} \right\|_2^2 + } \right.} }\\ {\left. {{\lambda _1}\left( {\sum\limits_{i = 1}^N {\sum\limits_{j = 1}^K {u_{ij}^m} } \left\| {\mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)} - {\mathit{\boldsymbol{c}}_j}} \right\|_2^2} \right)} \right) + }\\ {\frac{{{\lambda _2}}}{2}\sum\limits_{m = 1}^M {\left( {\left\| {{\mathit{\boldsymbol{W}}^{\left( m \right)}}} \right\|_{\rm{F}}^2 + \left\| {{\mathit{\boldsymbol{b}}^{\left( m \right)}}} \right\|_2^2} \right)} } \end{array} $

(14)

根据链式求导法则，分别对${\mathit{\boldsymbol{W}}^{(m)}}$和${\mathit{\boldsymbol{b}}^{(m)}}$求得梯度的表达式为

$ \frac{{\partial J}}{{\partial {\mathit{\boldsymbol{W}}^{\left( m \right)}}}} = \left( {{\Delta ^{\left( m \right)}} + {\lambda _1}{\Lambda ^{\left( m \right)}}} \right){\left( {{\mathit{\boldsymbol{h}}^{\left( {m - 1} \right)}}} \right)^{\rm{T}}} + {\lambda _2}{\mathit{\boldsymbol{W}}^{\left( m \right)}} $

(15)

$ \frac{{\partial J}}{{\partial {\mathit{\boldsymbol{b}}^{\left( m \right)}}}} = {\Delta ^{\left( m \right)}} + {\lambda _1}{\Lambda ^{\left( m \right)}} + {\lambda _2}{\mathit{\boldsymbol{b}}^{\left( m \right)}} $

(16)

由于重构误差和聚类误差同时作用在EDFCC算法的编码层，而解码层只用重构误差训练，因此可以将优化过程分为3种情况讨论，分别是:

1)$m = 1, \ldots, \frac{M}{2} - 1$

2)$m = \frac{M}{2}$

3)$m = \frac{M}{2} + 1, \ldots, M$

根据上述3种情况，将${\Delta ^{(m)}}$和${\Lambda ^{(m)}}$分别定义为

$ {\Delta ^{\left( m \right)}} = \left\{ \begin{array}{l} - \left( {{\mathit{\boldsymbol{x}}_i} - \mathit{\boldsymbol{h}}_i^{\left( M \right)}} \right) \odot g'\left( {\mathit{\boldsymbol{z}}_i^{\left( M \right)}} \right)\;\;\;\;\;\;m = M\\ {\left( {{\mathit{\boldsymbol{W}}^{\left( {m + 1} \right)}}} \right)^{\rm{T}}}{\Delta ^{\left( {m + 1} \right)}} \odot g'\left( {\mathit{\boldsymbol{z}}_i^{\left( M \right)}} \right)\;\;\;其他 \end{array} \right. $

(17)

$ \begin{array}{*{20}{c}} {{\Delta ^{\left( m \right)}} = }\\ {\left\{ \begin{array}{l} {\left( {{\mathit{\boldsymbol{W}}^{\left( {m + 1} \right)}}} \right)^{\rm{T}}}{\Delta ^{\left( {m + 1} \right)}} \odot g'\left( {\mathit{\boldsymbol{z}}_i^{\left( m \right)}} \right)\;\;\;\;\;\;\;\;\;\;m = 1, \cdots ,\frac{M}{2} - 1\\ \sum\limits_{j = 1}^k {u_{ij}^m} \left( {\mathit{\boldsymbol{h}}_i^{\left( {\frac{M}{2}} \right)} - {\mathit{\boldsymbol{c}}_j}} \right) \odot g'\left( {\mathit{\boldsymbol{z}}_i^{\left( M \right)}} \right)\;\;\;\;\;m = \frac{M}{2}\\ 0\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;m = \frac{M}{2} + 1, \cdots ,M \end{array} \right.} \end{array} $

(18)

式中，${\mathit{\boldsymbol{h}}_i}^{(0)} = {\mathit{\boldsymbol{x}}_i}$，符号⊙代表点乘，${\mathit{\boldsymbol{z}}_i}^{(m)}$和$g$′(·)分别表示为式(19)和式(20)

$ \mathit{\boldsymbol{z}}_i^{\left( m \right)} = {\mathit{\boldsymbol{W}}^{\left( m \right)}}\mathit{\boldsymbol{h}}_i^{\left( {m - 1} \right)} + {\mathit{\boldsymbol{b}}^{\left( m \right)}} $

(19)

$ \mathit{\boldsymbol{h'}}_i^{\left( m \right)} = g'\left( {{\mathit{\boldsymbol{W}}^{\left( m \right)}}\mathit{\boldsymbol{h}}_i^{\left( {m - 1} \right)} + {\mathit{\boldsymbol{b}}^{\left( m \right)}}} \right) $

(20)

使用随机梯度下降方法，DAE网络的权重$\mathit{\boldsymbol{W}}$和偏移量$\mathit{\boldsymbol{b}}$更新公式为

$ {\mathit{\boldsymbol{W}}^{\left( m \right)}} = {\mathit{\boldsymbol{W}}^{\left( m \right)}} - \mu \frac{{\partial J}}{{\partial {\mathit{\boldsymbol{W}}^{\left( m \right)}}}} $

(21)

$ {\mathit{\boldsymbol{b}}^{\left( m \right)}} = {\mathit{\boldsymbol{b}}^{\left( m \right)}} - \mu \frac{{\partial J}}{{\partial {\mathit{\boldsymbol{b}}^{\left( m \right)}}}} $

(22)

式中，参数$μ$代表学习率，用来控制模型收敛速度，一般设置为一个很小的值，过大时容易发生振荡现象。EDFCC算法的具体过程为:

输入：待聚类数据矩阵$\mathit{\boldsymbol{X}}\left({{\mathit{\boldsymbol{H}}^{\left(0 \right)}}} \right)$，聚类类别数$\mathit{\boldsymbol{K}}$。

输出：隶属度矩阵$\mathit{\boldsymbol{U}}$。

1) 初始化隶属度矩阵$\mathit{\boldsymbol{U}}$，并给出$\mathit{\boldsymbol{K}}$；

2) for $m$ = 1; $m$≤$M$; $m$++ do

3) 初始化${\mathit{\boldsymbol{W}}^{(m)}}$；

4) 初始化${\mathit{\boldsymbol{b}}^{(m)}}$；

5) end for;

6) while not convergence do

7) for $i$ = 1; $i$≤$M$; $i$++ do

8)${\mathit{\boldsymbol{h}}^{\left(i \right)}}$←(1)；

9) end for

10) 使用式(14)计算${J_{{\rm{rec}}}}$、${J_{\rm{f}}}$和${J_{{\rm{reg}}}}$；

11) 使用式(10)计算$\mathit{\boldsymbol{C}}$；

12) 对${\mathit{\boldsymbol{H}}^{\left({\frac{M}{2}} \right)}}$进行聚类并且使用式(9)更新$\mathit{\boldsymbol{U}}$；

13) for $m$ = 1; $m$≤$M$; $m$++ do

14) 使用式(21)更新${\mathit{\boldsymbol{W}}^{(m)}}$；

15) 使用式(22)更新${\mathit{\boldsymbol{b}}^{(m)}}$；

16) end for

17) end while

18) 返回隶属度矩阵$\mathit{\boldsymbol{U}}$。

为了验证EDFCC算法的有效性，使用Indian Pines和Pavia University数据集进行实验，并与k-means、SC、SSC、LRSC和AEKM 5种流行聚类算法的实验结果进行比较，采用准确率(ACC)和标准化互信息(NMI)等量化指标来验证结果的有效性，并且对EDFCC算法中$m$、${\lambda _1}$和激活函数的取值进行讨论。

美国印第安纳州农林数据集(Indian Pines)是由机载可视红外成像光谱仪(AVIRIS)获取的数据，其光谱范围为0.41~2.45 μm，空间分辨率为25 m，光谱分辨率为10 nm，共包含145×145个样本点，原始波段数有220个，除去水汽吸收带和低信噪比的波段，剩余200个波段作为研究的对象。该数据集共有16种地物类别，真实影像与真值图如图 3所示，具体的地物信息见表 1。

意大利帕维亚大学城区数据集(Pavia University)是由德国的机载反射光学光谱成像仪(ROSIS)获取的数据，其光谱范围为0.43~0.86 μm，空间分辨率为1.3 m，共包含610×340个样本点，原始波段数有115个，除去噪声波段，剩余103个波段作为研究的对象。该数据集共有9种地物类别，真实影像如图 4(a)所示，真值图如图 4(b)所示。考虑计算效率，选取该数据集的典型区域，共包含8种地物类别，典型区域的真值图如图 4(c)所示，具体的地物信息见表 2。

采用准确率(ACC)作为评价指标，衡量的是预测结果和真实标签是否一致，是一种简单而透明的评估手段，具体计算方法为

$ ACC\left( {\mathit{\boldsymbol{ \boldsymbol{\varOmega} }},\mathit{\boldsymbol{S}}} \right) = \frac{1}{N}\sum\limits_K {\mathop {\max }\limits_J } \left| {{\mathit{\boldsymbol{\omega }}_K} \cap {\mathit{\boldsymbol{S}}_J}} \right| $

(23)

式中，$\mathit{\boldsymbol{ \boldsymbol{\varOmega} }}=\left\{\boldsymbol{\omega}_{1}, \boldsymbol{\omega}_{2}, \cdots, \boldsymbol{\omega}_{K}\right\}$代表聚类结果的集合，$\boldsymbol{S}=\left\{\boldsymbol{S}_{1}, \boldsymbol{S}_{2}, \cdots, \boldsymbol{S}_{J}\right\}$代表高光谱图像原始分类的集合，$N$代表像元总数，$\boldsymbol{\omega}_{K}$代表预测的第$K$个簇，$\boldsymbol{S}_{J}$代表原始第$J$个类的集合。但是聚类结果仅能获得像元类别，不能确定此类别代表的含义以及类别间具有的联系或区别，并且ACC不能衡量聚类质量与聚类数量之间的关系，能衡量该关系的指标是标准化互信息(NMI)，它能够评估聚类结果与原始分类的关联程度，表示为

$ N\left( {\mathit{\boldsymbol{ \boldsymbol{\varOmega} }},\mathit{\boldsymbol{S}}} \right) = \frac{{I\left( {\mathit{\boldsymbol{ \boldsymbol{\varOmega} }},\mathit{\boldsymbol{S}}} \right)}}{{\left[ {H\left( \mathit{\boldsymbol{ \boldsymbol{\varOmega} }} \right) + H\left( \mathit{\boldsymbol{S}} \right)} \right]/2}} $

(24)

式中，$I$称为互信息，表示为

$ I\left( {\mathit{\boldsymbol{ \boldsymbol{\varOmega} }},\mathit{\boldsymbol{S}}} \right) = \sum\limits_K {\sum\limits_J {P\left( {{\mathit{\boldsymbol{\omega }}_k} \cap {\mathit{\boldsymbol{S}}_J}} \right){{\log }_2}\frac{{P\left( {{\mathit{\boldsymbol{\omega }}_k} \cap {\mathit{\boldsymbol{S}}_J}} \right)}}{{P\left( {{\mathit{\boldsymbol{\omega }}_k}} \right)P\left( {{\mathit{\boldsymbol{S}}_J}} \right)}}} } $

(25)

式中，$P({\mathit{\boldsymbol{\omega }}_K})$、$P({\mathit{\boldsymbol{S}}_J})$和$P({\mathit{\boldsymbol{\omega }}_K} \cap {\mathit{\boldsymbol{S}}_J})$分别表示像元属于聚类簇${\mathit{\boldsymbol{\omega }}_K}$的概率、属于原始分类${\mathit{\boldsymbol{S}}_J}$的概率以及属于聚类簇${\mathit{\boldsymbol{\omega }}_K}$和原始分类${\mathit{\boldsymbol{S}}_J}$交集的概率，H代表熵。

NMI的区间范围是0~1，0代表聚类结果与原始分类任意像元均不匹配，1代表聚类结果完全正确。NMI越大代表聚类效果越好。

实验的硬件环境为Intel i7-8550U处理器，8 GB内存；软件环境为Windows 10操作系统和MATLAB 2017b运算平台。为消除随机误差，在相同条件下进行了15次实验，取15次实验结果的聚类精度均值。为了比较不同参数的取值对实验结果带来的影响，固定${\lambda _2}$，以观察$m$、${\lambda _1}$和激活函数的不同取值对实验结果的影响。

本文算法与4种对比算法在两个数据集上的ACC和NMI见表 3和表 4。从表 3和表 4可以看出，相比于其他对比算法，本文提出的EDFCC算法具有最高的ACC值和NMI值，表明本文方法取得了最好的聚类精度。在Indian Pines数据集上，EDFCC算法的ACC和NMI相较于LRSC分别提升了0.03和0.01，在Pavia University数据集上EDFCC算法的ACC和NMI相较于LRSC分别提升了0.04和0.02，同时在两个数据集上相较于基于深度聚类的AEKM算法ACC分别提升了0.02和0.03。

k-means聚类算法具有最低的ACC值和NMI值，因为该算法不适用于处理非稀疏的高维数据，且只能发现球状簇；SC算法以及基于SC算法的SSC算法和LRSC算法的表现稍好，因为这些方法对数据分布的适应性更强，并且对数据进行了降维处理，因此在处理高维数据时比传统方法好。但是由于聚类过程中，需要选择欧氏距离、余弦相似度、高斯相似度等不同的相似性测度方式构建亲和矩阵，所以它们的聚类效果依赖于亲和矩阵，不同的亲和矩阵得到的聚类效果不同，并且无法提取数据的深层特征。

AEKM算法和本文提出的EDFCC算法都是基于深度学习的聚类算法，它们的效果优于传统聚类算法，这主要是因为深度神经网络具有降维能力、强大的非线性变换能力和深层特征表示学习的能力，使得基于深度学习的深度聚类算法模型的性能优于浅层聚类方法。同时本文提出的EDFCC算法获得了最好的聚类效果，优于AEKM算法，这是因为AEKM算法使用k-means聚类算法，而EDFCC算法使用FCM算法，软划分效果相较于硬划分稍好。

在效率上，本文提出的EDFCC算法运行时间介于k-means算法和SC、SSC、LRSC算法之间。由于EDFCC算法是基于深度学习的聚类算法，它在模型训练上要花费许多时间，所以相较于k-means算法效率要稍低；但是在降维效率上远远优于SC、SSC、LRSC算法，使得总体效率优于上述算法；由于软划分相较于硬划分效率稍低，使得与同样基于深度学习的AEKM算法相比总体效率稍低。

为了从直观上验证本文提出的EDFCC算法，以Pavia University数据集的典型区域为例，图 5给出了对比算法和EDFCC算法在其上的聚类效果。由图 5可以看出，k-means、SC、SSC和LRSC算法由于缺乏良好的非线性变化能力和深层特征表示学习能力，使得聚类时不同类别的像元混杂在一起，聚类结果图像元混杂程度高；AEKM和EDFCC算法通过应用非线性变换，使像元相对于新空间中相应的聚类中心变得紧凑，聚类结果图像元混杂程度大幅降低。EDFCC算法得到的聚类结果图噪声点最少，像元混杂程度最低。噪声点越少、像元混杂程度越低代表聚类越准确。这充分证明了本文提出的EDFCC算法的有效性。

系数${\lambda _2}$的目的是防止模型过拟合，固定${\lambda _2}$，观察${\lambda _1}$的不同取值对ACC和NMI的影响，从而得到最优的${\lambda _1}$值。如图 6所示，在大多数情况下，${\lambda _1}$的最优取值在0.1附近，说明此时模型重构过程与聚类过程平衡得最好。

FCM聚类算法是最常用的聚类算法之一，参数$m$控制着模糊类间的分享程度，所以使用FCM聚类算法时，选取恰当的权重指数$m$非常重要。虽然没有给出严格证明，但是在实际应用中，大多数研究者认为$m$的最佳取值范围是[1.5，2.5]，图 7展示了$m$在不同取值下ACC和NMI的变化，可以直观地看出当$m$在2.5附近时，ACC与NMI的值最高，这与经验相符，所以将$m$的最优值设定为2.5。

不同的激活函数导致模型的性能不同，使用的激活函数包括Tanh、Sigmoid、Relu和Softplus。EDFCC算法的ACC、NMI与激活函数的关系如图 8所示。由图 8可以看出，在实验过程中Tanh激活函数的效果优于其他3个激活函数，所以实验设置DAE网络的激活函数为Tanh。

为了有效利用高光谱图像的光谱信息，提出了一种嵌入式深度神经网络模糊C均值聚类算法(EDFCC)。该算法能够提取高光谱图像光谱信息的深层特征，提升聚类精度，促进高光谱图像聚类任务的研究。在高光谱图像聚类任务上，EDFCC算法利用FCM聚类算法约束特征映射过程，使得特征学习、参数优化和聚类同时进行，并且不需要额外的训练过程。实验表明，EDFCC算法明显优于其他聚类算法。然而基于DAE网络的聚类方法主要考虑重构原始输入，由于缺乏一定的高光谱图像空间信息约束，即使以聚类算法约束特征映射过程，聚类精度仍有待提升。下一步工作重点是对模型加以合适的高光谱图像空间信息约束，进一步提高模型的聚类精度。