自20世纪80年代以来，高光谱成像系统已成为具有广阔应用前景的智能医学辅助诊断工具，尤其是在疾病辅助诊断中(Akbari等，2010；Lu和Fei，2014；Lu和Mandal，2014；Carrión-Camacho等，2019)。高光谱图像(hyperspectral image, HSI)具有许多连续的窄光谱带(Fang等，2017；Wei等，2019；Zhang等，2021)，可以为区分生物组织的不同化学成分提供丰富的信息(Lu和Fei，2014)。然而，波段的增加导致数据包含大量冗余信息，对分类性能造成负面影响的同时增加了计算复杂度。降维(dimensionality reduction, DR)是解决此类问题的有效技术之一，该技术通过提取包含数据本征信息的低维特征以去除冗余信息。

研究者已提出许多可用于高光谱图像分析的降维方法。基于变换投影的降维方法根据特定的准则将原始空间特征投影到低维子空间中，以保证低维特征包含足够的有效信息。其中，经典的基于投影的线性方法包括主成分分析(principal component analysis, PCA)(Jolliffe和Cadima，2016)，Fisher线性判别分析(Fisher’ s linear discriminant analysis, LDA)(Yang和Yang，2003)和独立成分分析(independent component analysis, ICA)(Villa等，2011；Dalla等，2011)。为了有效挖掘数据固有的非线性结构，提出了局部保持投影(locality preserving projection, LPP)(He和Niyogi，2004)和邻域保持嵌入(neighborhood preserving embedding, NPE)(He等，2005)等。非参数加权特征提取(nonparametric weighted feature extraction, NWFE)(Kuo和Landgrebe，2004)、局部Fisher判别分析(local Fisher discriminant analysis, LFDA)(Li等，2012)、核LFDA (Li等，2011)和增强的Fisher判别分析(enhanced Fisher discriminant criterion, EFDC)(Gao等，2012)等，进一步将LPP与LDA集成在一起，以有效利用数据的局部信息和标签信息。

Yan等人(2007)基于图论提出了通用图嵌入(graph embedding, GE)框架，归纳了大多数的变换投影降维方法。在通用图嵌入框架中，无向图用于表征结构信息。目前，获得相对稀疏图的主要技术是$\varepsilon $球近邻和$k$近邻(k-nearest neighbor, kNN)。然而，这两种方法都对噪声敏感，并且在实践中通常难以设置合适参数。为了自适应地构造适当的图，引入稀疏表示(sparse representation, SR)(Wright等，2009, 2010)。通过探索数据的稀疏性，Cheng等人(2010)提出了一种稀疏图嵌入模型(sparse graph embedding model, SGE)。Ly等人(2014)利用类标签信息提出了基于稀疏图的判别分析(sparse graph-based discriminant analysis, SGDA)以提升SGE的性能。此外，为了降低计算复杂度，协同表示(collaborative representation, CR)(Zhang等，2011a)被用于构建替代图，并提出基于协作图的判别分析(collaborative graph-based discriminant analysis, CGDA)方法及CGDA的块形式(block version of CGDA, BCGDA)。通过挖掘同类和异类样本间的基于CR的数据流形结构，Lyu等人(2017)提出了判别流形图嵌入(collaborative discriminative manifold embedding, CDME)。

尽管SR和CR可以自适应地揭示数据局部结构，但它们忽略了全局结构。为了克服上述缺点，Chen等人(2012)在图构造过程中引入了低秩表示(low rank representation, LRR)。如果原始空间中的数据位于多个线性子空间的并集上，则LRR具有出色的性能(Sumarsono和Du，2015)。然而，现实中的数据并不总是线性的，特别是对于HSI数据，这可能导致LRR无法有效捕获几何结构(Cai等，2011)。为了保留局部几何结构，Liu等人(2014)提出了LapLRR(Laplacian regularized LRR)用于子空间聚类。另外，Li等人(2016)提出一种基于稀疏和低秩图的判别分析方法(sparse and low-rank graph based discriminant approach, SLGDA)及其核形式(Pan等，2017b)，以同时保留局部和全局结构信息。

然而，上述降维方法都是以向量为输入的，忽略了空间域中的重要相关性。尽管为克服此缺点研究了一些空间光谱联合分析技术(Tarabalka等，2010；Zhang等，2006, 2011b)，但仍将空间光谱特征转换为矢量进行分析，忽略了HSI的本质多维结构——张量。为了利用张量数据中的空间光谱信息，研究者已经采用多线性代数开发了基于张量的分析技术(Zhang等，2013；Bourennane等，2010；An等，2016)。局部张量判别分析通过引入张量块，在去除冗余信息的同时充分挖掘数据的空间信息(Zhang等，2013)。为了保留低维空间中的光谱和空间结构信息，Deng等人(2017, 2018)提出了张量局部保持投影(tensor LPP, TLPP)来学习投影矩阵。张量稀疏和低秩图的判别分析(tensor sparse and low-rank graph-based discriminant analysis, TSLGDA)通过在图构造过程中添加稀疏和低秩约束，获得了较好的降维结果(Pan等，2017a)。以上基于张量的方法已有效提升了分类精度，但它们缺乏对其他结构特征的有效挖掘(例如固有几何信息)。

在张量分析中，增强判别能力和利用固有结构信息是两个重要的问题。但是，现有的基于张量表示的降维方法通常分开考虑这两个问题，并且仅利用部分的固有结构信息。本文通过构建一个集成的框架解决这些问题。综合考虑张量分析，多流形信息和判别信息，设计了一种基于张量的拉普拉斯稀疏低秩图(tensor-based Laplacian regularized sparse and low-rank graph, T-LapSLRG)嵌入技术，用于医学HSI的判别分析。在T-LapSLRG中，利用张量表示来捕获空间信息，同时引入稀疏和低秩约束以利用数据的局部和全局结构。多流形正则化约束用以学习具有判别信息的低维特征，并进一步调整基于张量的权重图以防止目标过拟合，这使得所提出的T-LapSLRG比其他基于张量的降维方法更具鲁棒性和鉴别力。本文首次将基于张量图的多流形图嵌入技术应用于医学高光谱图像分析。

慢性肾脏病(chronic kidney disease，CKD)是全球性的公共卫生问题(徐金升等，2018)，其发病率超过10 %。在CKD中，膜性肾病(membranous nephropathy, MN)是导致成人肾病综合征的一个常见病因(Ronco和Debiec，2015；van den Brand等，2011；徐金升等，2018)。据报道，肾活检标本中MN的比例在全世界范围内持续升高(Dong等，2016；Schweitzer等，2015)。膜性肾病的传统诊断方法主要依靠血清学特征和肾脏病理学特征，已经不能满足临床对准确性和智能性的诊断需求。此外，传统光学图像上原发膜性肾病(primary MN，PMN)和乙肝相关膜性肾病(hepatitis B virus-associated MN，HBV-MN)之间的病理特征高度相似，很难用光学显微镜区分。因此本文验证数据采用中日友好医院采集的膜性肾炎数据，其中包含肾脏病科的19位不同患者的肾活检组织切片的54幅显微高光谱图像(10名HBV-MN患者的30幅图像和9名PMN患者的24幅图像)。实验结果表明，所提出的T-LapSLRG在膜性肾病诊断中具有潜在的临床应用价值。

张量用于定义多维数组(Zhang等，2013)。多维数组的维数称为张量的阶数。对于$N$阶张量${\boldsymbol{\mathcal{X}}} \in {{\bf{R}}^{{I_1} \times \ldots \times {I_n} \times \ldots \times {I_N}}}$，其元素表示为${{\boldsymbol{\mathcal{X}}}_{{i_1} \ldots {i_n} \ldots {i_N}}}$，其中$1{\rm{ }} \le {i_n} \le {I_n}(1{\rm{ }} \le n \le N)$。下面介绍几个相关的张量基本原理和定义。

定义1   $k$模展开。假设${\boldsymbol{\mathcal{X}}} \in {{\bf{R}}^{{I_1} \times \ldots \times {I_n} \times \ldots \times {I_N}}}$是一个$N$阶张量，${{\boldsymbol{X}}^{(k)}} \in {\bf{R}}{I_k} \times ({I_1} \ldots {I_{k - 1}}{I_{k + 1}} \ldots {I_N})$是$\boldsymbol{\mathcal{X}}$的$k$模展开，它由$\boldsymbol{\mathcal{X}}$的所有$k$模向量按列组成。

定义2  $k$模积。将${\boldsymbol{\mathcal{Y}}} = {\boldsymbol{\mathcal{X}}}{ \times _k}{\bf{U}}$表示为张量$\boldsymbol{\mathcal{X}}$和矩阵${\boldsymbol{U}} \in {{\boldsymbol{R}}^{J \times {I_k}}}$的$k$模积，其中${\boldsymbol{\mathcal{Y}}} = {{\bf{R}}^{{I_1} \times \cdots \times {I_{k - 1}} \times J \times {I_{k + 1}} \times \cdots \times {I_N}}}$，$\boldsymbol{\mathcal{Y}}$中每一项计算为

$ \boldsymbol{\mathcal{Y}}_{i_{1} \cdots i_{k-j} ji_{k+1} \cdots i_{N}}=\sum\limits_{i_{k}=1}^{I_{k}} \boldsymbol{\mathcal{X}}_{i_{1} i_{2} \cdots i_{N}} \boldsymbol{U}_{ji_{k}} $

(1)

张量的展开可用于表示$k$模积，即

$ \boldsymbol{\mathcal{Y}}=\boldsymbol{\mathcal{X}} \times{ }_{k} \boldsymbol{U} \Leftrightarrow \boldsymbol{\mathcal{Y}}^{(k)}=\boldsymbol{U} \boldsymbol{\mathcal{X}}^{(k)} $

(2)

定义3  F-范数(Frobenius范数)。张量$\boldsymbol{\mathcal{X}}$的F-范数为$\left| {\boldsymbol{\mathcal{X}}} \right|_{\rm{F}} = {\left({\sum\limits_{{i_1} \cdots {i_N}} {{{\left({{{\boldsymbol{\mathcal{X}}}_{{i_1} \cdots {i_N}}}} \right)}^2}} } \right)^{1/2}}$。

定义4  张量缩并。张量${\boldsymbol{\mathcal{X}}} \in {{\bf{R}}^{{I_1} \times \cdots \times {I_N} \times f{'_1} \times \cdots \times I{'_{N'}}}}$和张量${\boldsymbol{\mathcal{Y}}} \in {{\bf{R}}^{{I_1} \times \cdots \times {I_N} \times \cdots \times I'{'_1} \times \cdots \times I'{'_{N''}}}}$的缩并定义为

$ \begin{aligned} &{[\boldsymbol{\mathcal{X}} \otimes \boldsymbol{\mathcal{Y}} ;(1: N)(1: N)]_{i_{1}, i_{2}, \cdots, i_{N}}=} \\ &\sum\limits_{i_{1}=1}^{I_{1}} \cdots \sum\limits_{i_{N}=1}^{I_{N}} \boldsymbol{\mathcal{X}}_{i_{1}, \cdots, i_{N}, i_1^{\prime}, \cdots, i_{N^{\prime}}^{\prime}} \boldsymbol{\mathcal{Y}}_{i_{1}, \cdots, i_{N}, i_1^{\prime \prime}, \cdots, i_{N^{\prime \prime}}^{\prime \prime}} \end{aligned} $

(3)

式中，$[{\boldsymbol{\mathcal{X}}} \otimes {\boldsymbol{\mathcal{Y}}};(1:N)(1:N)] \in {{\bf{R}}^{I{'_1} \times \cdots \times I{'_{N'}} \times I'{'_1} \times \cdots \times I'{'_{N''}}}}$。张量缩并的必要条件是两个张量维度相同。特别地，当张量缩并中张量$\boldsymbol{\mathcal{X}}$, ${\boldsymbol{\mathcal{Y}}} \in {{\bf{R}}^{{I_1} \times \cdots \times {I_k} \times \cdots \times {I_N}}}$维度的所有下标中不包含$k$时，该缩并可表示为$[{\boldsymbol{\mathcal{X}}} \otimes {\boldsymbol{\mathcal{Y}}};(\bar k)(\bar k)]$，且有

$ [\boldsymbol{\mathcal{X}} \otimes \boldsymbol{\mathcal{Y}} ;(\bar{k})(\bar{k})]=\boldsymbol{\mathcal{X}}^{(k)} \boldsymbol{\mathcal{Y}}^{(k)^{\mathrm{T}}} $

(4)

医学高光谱图像可自然地视为一个三阶张量${\boldsymbol{\mathcal{T}}} \in {{\bf{R}}^{{I_1} \times {I_2} \times {I_3}}}$，其中, ${I_1}$和${I_2}$表示立方体的空间维尺度，${I_3} = d$表示光谱维度。从其中抽取的第$k$个张量块可表示为${{\boldsymbol{\mathcal{X}}}_k} \in {{\bf{R}}^{{i_1} \times {i_2} \times d}}$，由第$k$个像素及其所在空间域中的${i_1} \times {\rm{ }}{i_2}$个邻接像素组成。假设医学高光谱图像的感兴趣区域包含$M$个像素($M \le {I_1} \times {I_2}$)，则由$M$个标记块构成的训练集$\left\{ {{{\boldsymbol{\mathcal{X}}}_k}} \right\}_{k = 1}^M$可表示为四阶张量${{\boldsymbol{\mathcal{X}}}^{(*)}} \in {{\bf{R}}^{{i_1} \times {i_2} \times d \times M}}$。第$c$类的张量块集合可表示为$\left\{ {{{\boldsymbol{\mathcal{X}}}_{k, c}}} \right\}_{k = 1}^{{M_c}}$，也可以表示为四阶张量${{\boldsymbol{\mathcal{X}}}^{(c)}} \in {{\bf{R}}^{{i_1} \times {i_2} \times d \times {M_c}}}$的形式，其中${M_c}$表示第$c$类包含张量块的个数，且满足$\sum\limits_{c = 1}^C {{M_c}} = M(c \in \{ 1, 2, \cdots, C\})$。

为了充分利用数据的标签信息、固有几何结构、低秩和稀疏特性进行高光谱图像分析，构建基于张量的拉普拉斯稀疏低秩图模型

$ \begin{gathered} \min \limits_{W(c)} \frac{1}{2}\left\|\boldsymbol{\mathcal{X}}^{(c)}-\boldsymbol{\mathcal{X}}^{(c)} \times{ }_{4} \boldsymbol{W}^{(c)}\right\|_{\mathrm{F}}^{2}+\lambda\left\|\boldsymbol{W}^{(c)}\right\|_{*}+ \\ \beta\left\|\boldsymbol{W}^{(c)}\right\|_{1}+\frac{\gamma}{2} \operatorname{tr}\left(\boldsymbol{W}^{(c) \mathrm{T}} \boldsymbol{Z}^{(c)} \boldsymbol{W}^{(c)}\right) \\ \text { s. t. } \operatorname{diag}\left(\boldsymbol{W}^{(c)}\right)=0 \end{gathered} $

(5)

式中, $\lambda, \beta $和$\gamma $是用于控制各项比重的平衡参数，${{\boldsymbol{W}}^{(c)}} \in {{\bf{R}}^{{M_c} \times {M_c}}}$表示第$c$类张量块的权重矩阵，${{\boldsymbol{Z}}^{(c)}}$表示仿射矩阵${{\boldsymbol{A}}^{(c)}}$的拉普拉斯矩阵, ${ \times _4}$运算是4模积(定义2)。矩阵${{\boldsymbol{A}}^{(c)}}$的第$p$行，第$q$列对应元素的计算式为

$ \boldsymbol{A}_{p q}^{(c)}=\exp \left(-\left\|{\boldsymbol{\mathcal{X}}}_{p, c}-\boldsymbol{\mathcal{X}}_{q, c}\right\|_{\mathrm{F}}^{2} / \eta_{p} \eta_{q}\right) $

(6)

式中, ${\eta _p} = {\left\| {{{\boldsymbol{\mathcal{X}}}_{p, c}} - {\cal X}_{p, c}^{(kn)}} \right\|_{\rm{F}}}$是第$c$类中第$p$个样本${{\boldsymbol{\mathcal{X}}}_{p, c}}$的近邻元素的尺度，${\boldsymbol{\mathcal{X}}}_{p, c}^{({\rm{knn}})}$表示${\boldsymbol{\mathcal{X}}}_{p, c}$的第$k$个最近邻。全局权重矩阵${\boldsymbol{W}}$可表示为

$ \boldsymbol{W}=\left[\begin{array}{ccc} \boldsymbol{W}^{(1)} & \cdots & 0 \\ \vdots & \ddots & \vdots \\ 0 & \cdots & \boldsymbol{W}^{(C)} \end{array}\right] $

(7)

采用交替方向乘子法(alternating direction method of multipliers, ADMM)解决优化问题(5)。为了使目标函数可分离，引入两个辅助变量${{\boldsymbol{Q}}^{(c)}}$, ${{\boldsymbol{J}}^{(c)}}$。优化问题(5)可转化为

$ \begin{gathered} \min \limits_{\boldsymbol{W}^{(c)}} \frac{1}{2}\left\|\boldsymbol{\mathcal{X}}^{(c)}-\boldsymbol{\mathcal{X}}^{(c)} \times{ }_{4} \boldsymbol{W}^{(c)}\right\|_{\mathrm{F}}^{2}+\lambda\left\|\boldsymbol{Q}^{(c)}\right\|_{*}+ \\ \qquad \beta\left\|\boldsymbol{J}^{(c)}\right\|_{1}+\frac{\gamma}{2} \operatorname{tr}\left(\boldsymbol{W}^{(c) \mathrm{T}} \boldsymbol{Z}^{(c)} \boldsymbol{W}^{(c)}\right) \\ \text { s. t. } \quad \boldsymbol{W}^{(c)}=\boldsymbol{Q}^{(c)}, \boldsymbol{W}^{(c)}=\boldsymbol{J}^{(c)}-\operatorname{diag}\left(\boldsymbol{J}^{(c)}\right) \end{gathered} $

(8)

则增广拉格朗日函数为

$ \begin{gathered} L\left(\boldsymbol{Q}^{(c)}, \boldsymbol{J}^{(c)}, \boldsymbol{W}^{(c)}, \boldsymbol{D}_{1}, \boldsymbol{D}_{2}\right)= \\ \frac{1}{2}\left\|\boldsymbol{\mathcal{X}}^{(c)}-\boldsymbol{\mathcal{X}}^{(c)} \times{ }_{4} \boldsymbol{W}^{(c)}\right\|_{\mathrm{F}}^{2}+\lambda\left\|\boldsymbol{Q}^{(c)}\right\|_{*}+\beta\left\|\boldsymbol{J}^{(c)}\right\|_{1}+ \\ \frac{\gamma}{2} \operatorname{tr}\left(\boldsymbol{W}^{(c) \mathrm{T}} \boldsymbol{Z}^{(c)} \boldsymbol{W}^{(c)}\right)+\left\langle\boldsymbol{D}_{1}, \boldsymbol{W}^{(c)}-\boldsymbol{Q}^{(c)}\right\rangle+ \\ \left\langle\boldsymbol{D}_{2}, \boldsymbol{W}^{(c)}-\boldsymbol{J}^{(c)}+\operatorname{diag}\left(\boldsymbol{W}^{(c)}\right)\right\rangle+ \\ \frac{\mu}{2}\left(\left\|\boldsymbol{W}^{(c)}-\boldsymbol{Q}^{(c)}\right\|_{\mathrm{F}}^{2}+\left\|\boldsymbol{W}^{(c)}-\boldsymbol{J}^{(c)}+\operatorname{diag}\left(\boldsymbol{J}^{(c)}\right)\right\|_{\mathrm{F}}^{2}\right) \end{gathered} $

(9)

式中，$\mu $为参数，${{\boldsymbol{D}}_1}$, ${{\boldsymbol{D}}_2}$为拉格朗日算子。

为了使函数$L\left({{{\boldsymbol{Q}}^{(c)}}, {{\boldsymbol{J}}^{(c)}}, {{\boldsymbol{W}}^{(c)}}, {{\boldsymbol{D}}_1}, {{\boldsymbol{D}}_2}} \right)$最小化，通过固定剩余变量来交替更新每个变量。更新规则以如下方法实现

$ \begin{gathered} \boldsymbol{Q}_{t+1}^{(c)}=\arg \underset{\boldsymbol{Q}^{(c)}}{\min } \lambda\left\|\boldsymbol{Q}^{(c)}\right\|_{*}+\left\langle\boldsymbol{D}_{1, t}, \boldsymbol{W}_{t}^{(c)}-\boldsymbol{Q}^{(c)}\right\rangle+\\ \frac{\mu_{t}}{2}\left\|\boldsymbol{W}_{t}^{(c)}-\boldsymbol{Q}^{(c)}\right\|_{\mathrm{F}}^{2}=\arg \min \limits_{Q^{(c)}} \frac{\lambda}{\mu_{t}}\left\|\boldsymbol{Q}^{(c)}\right\|_{*}+ \\ \frac{1}{2}\left\|\boldsymbol{Q}^{(c)}-\left(\boldsymbol{W}_{t}^{(c)}+\frac{\boldsymbol{D}_{1, t}}{\mu_{t}}\right)\right\|_{\mathrm{F}}^{2}=\boldsymbol{\varOmega}_{\frac{\lambda}{\mu_{t}}}\left(\boldsymbol{W}_{t}^{(c)}+\frac{\boldsymbol{D}_{1, t}}{\mu_{t}}\right) \end{gathered} $

(10)

$ \begin{gathered} \boldsymbol{J}_{t+1}^{(c)}=\arg \min \limits_{\boldsymbol{J}^{(c)}} \beta\left\|\boldsymbol{J}^{(c)}\right\|_{1}+\left\langle\boldsymbol{D}_{2, t}, \boldsymbol{W}_{t}^{(c)}-\boldsymbol{J}^{(c)}\right\rangle+ \\ \frac{\mu_{t}}{2}\left\|\boldsymbol{W}_{t}^{(c)}-\boldsymbol{J}^{(c)}\right\|_{\mathrm{F}}^{2}=\arg \min \limits_{\boldsymbol{J}^{(c)}} \frac{\beta}{\mu_{t}}\left\|\boldsymbol{J}^{(c)}\right\|_{1}+ \\ \frac{1}{2}\left\|\boldsymbol{J}^{(c)}-\left(\boldsymbol{W}_{t}^{(c)}+\frac{\boldsymbol{D}_{2, t}}{\mu_{t}}\right)\right\|_{\mathrm{F}}^{2}=\boldsymbol{\mathcal{S}}_{\frac{\beta}{\mu_{t}}}\left(\boldsymbol{W}_{t}^{(c)}+\frac{\boldsymbol{D}_{2, t}}{\mu_{t}}\right), \\ \boldsymbol{J}_{t+1}^{(c)}=\boldsymbol{J}_{t+1}^{(c)}-\operatorname{diag}\left(\boldsymbol{J}_{t+1}^{(c)}\right) \end{gathered} $

(11)

式中, ${\mu _t}$表示学习率，${\mathit{\Omega }_\tau }(\Delta) = {\boldsymbol{Q}}{{\boldsymbol{\mathcal{S}}}_\tau }(\Sigma){{\boldsymbol{V}}^{\rm{T}}}$表示奇异值阈值算子，并且软阈值算子表示为${{\boldsymbol{\mathcal{S}}}_\tau }(x) = {\mathop{\rm sgn}} (x)\max (|x| - \tau, 0)$(Cai等，2010)。${\boldsymbol{W}}_{t + 1}^{(c)}$可以在${\boldsymbol{Q}}_{t + 1}^{(c)}{\rm{ }}$和${\boldsymbol{J}}_{t + 1}^{(c)}$固定的情况下获得，即

$ \begin{gathered} \boldsymbol{W}_{t+1}^{(c)}=\arg \min \limits_{\boldsymbol{W}^{(c)}} \frac{1}{2}\left\|{\boldsymbol{\mathcal{X}}}^{(c)}-{\boldsymbol{\mathcal{X}}}^{(c)} \times{ }_{4} \boldsymbol{W}^{(c)}\right\|_{\mathrm{F}}^{2}+ \\ \frac{\gamma}{2} \operatorname{tr}\left(\boldsymbol{W}^{(c)^{\mathrm{T}}} \boldsymbol{Z}^{(c)} \boldsymbol{W}^{(c)}\right)+\left\langle\boldsymbol{D}_{1, t}, \boldsymbol{W}^{(c)}-\boldsymbol{Q}_{t+1}^{(c)}\right\rangle+ \\ \left\langle\boldsymbol{D}_{2, t}, \boldsymbol{W}^{(c)}-\boldsymbol{J}_{t+1}^{(c)}\right\rangle+\frac{\mu_{t}}{2}\left(\left\|\boldsymbol{W}^{(c)}-\boldsymbol{Q}_{t+1}^{(c)}\right\|_{\mathrm{F}}^{2}+\right. \\ \left.\left\|\boldsymbol{W}^{(c)}-\boldsymbol{J}_{t+1}^{(c)}\right\|_{\mathrm{F}}^{2}\right)=\left(\boldsymbol{M}^{(c)}+2 \mu_{t} \boldsymbol{I}+\gamma \boldsymbol{Z}^{(c)}\right)^{-1} \times \\ \left(\boldsymbol{M}^{(c)}+\mu_{t} \boldsymbol{Q}_{t+1}^{(c)}+\mu_{t} \boldsymbol{J}_{t+1}^{(c)}-\boldsymbol{D}_{1, t}-\boldsymbol{D}_{2, t}\right) \end{gathered} $

(12)

式中，${\boldsymbol{I}} \in {{\bf{R}}^{{M_c} \times {M_c}}}$是单位矩阵，并且${{\boldsymbol{M}}^{(c)}} = \left[ {{{\boldsymbol{\mathcal{X}}}^{(c)}} \otimes } \right.\left. {{{\boldsymbol{\mathcal{X}}}^{(c)}};(\bar 4)(\bar 4)} \right]$。

在计算所有子相似度矩阵后，构造全局权重矩阵${\boldsymbol{W}}$。利用集合和权重矩阵${\boldsymbol{W}}$构造一种新型的基于张量的拉普拉斯稀疏低秩图${\boldsymbol{G}} = \{ {\boldsymbol{\mathcal{X}}}, {\boldsymbol{W}}\} $。

基于张量的图嵌入(tensor-based graph embeding, T-GE)的降维方法是找到变换矩阵$\{ {{\boldsymbol{U}}_1}, {{\boldsymbol{U}}_2}, \ldots, {{\boldsymbol{U}}_N}\} $，利用变换矩阵将高维样本${\boldsymbol{\mathcal{X}}_i}$投影到低维子空间$\tilde{\boldsymbol{\mathcal{X}}}_{i}$中，其中${\tilde{\boldsymbol{\mathcal{X}}}_{i}} = {{\boldsymbol{\mathcal{X}}}_i}{ \times _1}{{\boldsymbol{U}}_1}{ \times _2}{{\boldsymbol{U}}_2} \cdots { \times _N}{{\boldsymbol{U}}_N}$。目的是在高维空间中保留图${\boldsymbol{G}} = \{ {\boldsymbol{\mathcal{X}}}, {\boldsymbol{W}}\} $包含的结构信息。

因此，将T-GE的优化问题构造为

$ \begin{gathered} \min \limits_{\boldsymbol{U}_{1}, \boldsymbol{U}_{2}, \cdots, \boldsymbol{U}_{N}} \sum\limits_{i, j}\left\|\tilde{\boldsymbol{\mathcal{X}}}_{i}-\tilde{\boldsymbol{\mathcal{X}}}_{j}\right\|^{2} \boldsymbol{W}_{i j} \\ \text { s. t. } \sum\limits_{i, j}\left\|\boldsymbol{\mathcal{X}}_{i} \times_{1} \boldsymbol{U}_{1} \times{ }_{2} \boldsymbol{U}_{2} \cdots \times_{N} \boldsymbol{U}_{N}\right\|^{2} \boldsymbol{B}_{i i}=1 \end{gathered} $

(13)

式中, ${{\boldsymbol{B}}_{ii}} = \sum\limits_j {{{\boldsymbol{W}}_{ij}}} $。

上述优化问题采用迭代方案来解决。假设$\{ {{\boldsymbol{U}}_1}, \ldots, {{\boldsymbol{U}}_{k - 1}}, {{\boldsymbol{U}}_{k + 1}}, \ldots, {{\boldsymbol{U}}_N}\} $是固定的，令${\boldsymbol{\mathcal{X}}_{i, (k)}} = {\boldsymbol{\mathcal{X}}_i}{ \times _1}{{\boldsymbol{U}}_1} \ldots { \times _{k - 1}}{{\boldsymbol{U}}_{k - 1}}{ \times _{k + 1}}{{\boldsymbol{U}}_{k + 1}} \ldots {\rm{ }}{ \times _N}{{\boldsymbol{U}}_N}$，优化问题(13)转化为

$ \begin{gathered} \min \limits_{\boldsymbol{U}_{k}} J\left(\boldsymbol{U}_{k}\right)= \\ \sum\limits_{i, j}\left\|\tilde{\boldsymbol{\mathcal{X}}}_{i, (k)} \times{ }_{k} \boldsymbol{U}_{k}-\tilde{\boldsymbol{\mathcal{X}}}_{j, (k)} \times{ }_{k} \boldsymbol{U}_{k}\right\|^{2} \boldsymbol{W}_{i j}= \\ \sum\limits_{i, j}\left\|\boldsymbol{U}_{k} \tilde{\boldsymbol{X}}_{i}^{k}-\boldsymbol{U}_{k} \tilde{\boldsymbol{X}}_{j}^{k}\right\|^{2} \boldsymbol{W}_{i j}= \\ \sum\limits_{i, j} \operatorname{tr}\left(\boldsymbol{U}_{k}\left(\left(\tilde{\boldsymbol{X}}_{i}^{k}-\tilde{\boldsymbol{X}}_{j}^{k}\right)\left(\tilde{\boldsymbol{X}}_{i}^{k}-\tilde{\boldsymbol{X}}_{j}^{k}\right)^{\mathrm{T}} \boldsymbol{W}_{i j}\right) \boldsymbol{U}_{k}^{\mathrm{T}}\right)= \\ \operatorname{tr}\left(\boldsymbol{U}_{k}\left(\sum\limits_{i, j}\left(\tilde{\boldsymbol{X}}_{i}^{k}-\tilde{\boldsymbol{X}}_{j}^{k}\right)\left(\tilde{\boldsymbol{X}}_{i}^{k}-\tilde{\boldsymbol{X}}_{j}^{k}\right)^{\mathrm{T}} \boldsymbol{W}_{i j}\right) \boldsymbol{U}_{k}^{\mathrm{T}}\right) \\ \text { s. t. } \operatorname{tr}\left(\boldsymbol{U}_{k}\left(\sum\limits_{i} \tilde{\boldsymbol{X}}_{i}^{k}\left(\tilde{\boldsymbol{X}}_{i}^{k}\right)^{\mathrm{T}} \boldsymbol{B}_{i i}\right) \boldsymbol{U}_{k}^{\mathrm{T}}\right)=1 \end{gathered} $

(14)

式中, $\tilde{\boldsymbol{X}}_{i}^{k}$是$\tilde{\boldsymbol{\mathcal{X}}}_{i}^{k}$的$k$模展开。

通过广义特征值问题

$ \begin{gathered} \left(\sum\limits_{i, j}\left(\tilde{\boldsymbol{X}}_{i}^{k}-\tilde{\boldsymbol{X}}_{j}^{k}\right)\left(\tilde{\boldsymbol{X}}_{i}^{k}-\tilde{\boldsymbol{X}}_{j}^{k}\right)^{\mathrm{T}} \boldsymbol{W}_{i j}\right) \boldsymbol{u}= \\ \lambda\left(\sum\limits_{i} \tilde{\boldsymbol{X}}_{i}^{k}\left(\tilde{\boldsymbol{X}}_{i}^{k}\right)^{\mathrm{T}} \boldsymbol{B}_{i i}\right) \boldsymbol{u} \end{gathered} $

(15)

获取的最小$k$个特征值对应的特征向量构成了最佳投影矩阵${{\boldsymbol{U}}_k}$。算法1中给出了T-LapSLRG的详细描述，算法中${\mu _0}$为参数$\mu $的初始值，${\mu _{{\rm{max}}}}$为参数$\mu $的最大值，${\rho _0}$为权重常数，${\varepsilon _1}$和${\varepsilon _2}$是较小的阈值常数，t为迭代次数。

算法1  基于张量的拉普拉斯稀疏低秩图嵌入

输入：有标签训练数据$\left\{ {{{\boldsymbol{\mathcal{X}}}_k}} \right\}_{k = 1}^M$，参数$\lambda, \beta $和$\gamma $，张量块的空间尺寸${i_1} \times {\rm{ }}{i_2}$和降低维度集合$\{ {B_1}, {B_2}, \ldots, {B_N}\} $。

初始化：${\boldsymbol{Q}}_0^{(c)} = {\boldsymbol{J}}_0^{(c)} = {\boldsymbol{W}}_0^{(c)} = {\bf{0}}$, ${{\boldsymbol{D}}_{1, 0}} = {{\boldsymbol{D}}_{2, 0}} = {\bf{0}}$, ${\mu _0} = 0.1$, ${\mu _{\max }} = {10^3}$, ${\rho _0} = 1.1$, ${\varepsilon _1} = {10^{ - 4}}$, ${\varepsilon _2} = {10^{ - 3}}$, $t = 0$, $\mathit{maxIter} = 100$。

从$c{\rm{ }} = 1, 2, \ldots, C$开始循环，根据式(10)—式(12)，计算${\boldsymbol{Q}}_{t + 1}^{(c)}$, ${\boldsymbol{J}}_{t + 1}^{(c)}$和${\boldsymbol{W}}_{t + 1}^{(c)}$。

更新拉格朗日乘子：

$ \boldsymbol{D}_{1, t+1}=\boldsymbol{D}_{1, t}+\mu_{t}\left(\boldsymbol{W}_{t+1}^{(c)}-\boldsymbol{Q}_{t+1}^{(c)}\right), $

$ \boldsymbol{D}_{2, t+1}=\boldsymbol{D}_{2, t}+\mu_{t}\left(\boldsymbol{W}_{t+1}^{(c)}-\boldsymbol{J}_{t+1}^{(c)}\right)。$

更新$\mu :{\mu _{t + 1}} = \min \left({\rho {\mu _t}, {\mu _{\max }}} \right)$。式中，

$ \rho= \begin{cases}\rho_{0} & \frac{\mu_{t} \max \left(\begin{array}{l} \left\|\boldsymbol{W}_{t+1}^{c}-\boldsymbol{W}_{t}^{c}\right\|_{\mathrm{F}}, \\ \left\|\boldsymbol{Q}_{t+1}^{c}-\boldsymbol{Q}_{t}^{c}\right\|_{\mathrm{F}}, \\ \left\|\boldsymbol{J}_{t+1}^{c}-\boldsymbol{J}_{t}^{c}\right\|_{\mathrm{F}} \end{array}\right)}{\left\|\boldsymbol{\mathcal{X}}^{c}\right\|_{\mathrm{F}}}<\varepsilon_{2} \\ 0 & {\text { 其他 }}\end{cases} $

检查收敛条件

$ \begin{aligned} &\left\|\boldsymbol{W}_{t+1}^{(c)}-\boldsymbol{Q}_{t+1}^{(c)}\right\|_{\infty}<\varepsilon_{1}, \left\|\boldsymbol{W}_{t+1}^{(c)}-\boldsymbol{J}_{t+1}^{(c)}\right\|_{\infty}<\varepsilon_{1} \\ &t \longleftarrow t+1 \end{aligned} $

当满足收敛条件或$t > \mathit{maxIter}$结束循环。

根据式(7)构造矩阵${\boldsymbol{W}}$。

根据式(15)计算投影矩阵${{\boldsymbol{U}}_k}, k = 1, \cdots, N$。

输出：投影矩阵集$\{ {{\boldsymbol{U}}_1}, {{\boldsymbol{U}}_2}, \ldots, {{\boldsymbol{U}}_N}\} $。

利用膜性肾病医学高光谱数据集来评估所提出技术的有效性。为了捕获肾活检组织切片的显微高光谱图像，将推扫式高光谱成像系统SOC-710与生物显微镜CX31RTSF结合构造显微高光谱成像系统(见图 1)。SOC-710系统捕获具有空间尺寸696×520的128个光谱带，以4.69 nm的分辨率覆盖了400~1 000 nm的光谱范围。

该数据集包含从中日友好医院肾脏病科的19位不同患者的肾活检组织切片拍摄的54幅显微高光谱图像。肾脏病理切片用过碘酸六胺银+马松染色(PASM + M)染色，以准确显示和定位肾小球疾病中的各种免疫复合物。获得的医学高光谱图像(medical HSI, MHSI)数据集由10名HBV-MN患者的30幅HBV-MN图像和9名PMN患者的24幅PMN图像组成。为了获得足够的信息，从每位患者处收集了3幅图像。但是，基于每个图像都包含完整的肾小球的图像获取原理，其中2位PMN患者的图像少于3幅。经验丰富的专科肾脏病专家标记了图像上免疫复合物的区域。图 2展示了分别从HBV-MN患者和PMN患者获得的伪彩色图像和标记图像的示例, 其中白色区域为感兴趣物质标注区域。验证数据的详细信息在表 1中列出。加粗标记的患者是由临床医生诊断为具有典型临床症状的患者。因此，将PMN患者(ID：15684)和HBV-MN患者(ID：17472)用作后续实验的训练数据。

T-LapSLRG方法包含5个参数，包括窗口大小${i_1} \times {\rm{ }}{i_2}$，平衡参数$\lambda $，$\beta $和$\gamma $以及降低维度集$\{ {B_1}, {B_2}, \cdots, {B_N}\} $。$\lambda $，$\beta $和$\gamma $分别控制低秩项、稀疏项和流形正则项，参数范围为{10^-4, 10^-3, 10^-2, 10^-1, 10⁰, 10¹, 10²}。为了评估T-LapSLRG的有效性，引入了经典支持向量机分类器(support vector machine, SVM)(Moughal，2013)对提取的低维特征进行分类。分类性能通过总体准确性(overall accuracy, OA)进行评估。经初步实验验证，$\lambda $和$\beta $的影响是相当的，因此展示了$\lambda = \beta $时的实验结果。图 3说明了MN数据集上具有不同$\lambda $，$\beta $和$\gamma $值的OA。结果表明，将$\lambda $, $\beta $和$\gamma $设置为(10^-1, 10^-1, 10¹)可获得最佳性能。因此，在以下实验中，选择(10^-1, 10^-1, 10¹)作为权衡参数。

假设像素位于窗口的中心，使用固定的空间窗口来提取与每个像素相对应的张量样本，即${i_1} = {i_2} = w$。每个张量的标签视为其对应中心像素的标签。窗口的大小起着非常重要的作用，随着窗口大小的增加，张量样本中包含的像素将增加，这将提供更多的空间光谱信息。但是，过大的窗口可能会导致张量样本包含冗余信息。实际上，窗口的尺寸通常根据经验和实验确定。表 2列出了以(3×3, 5×5, 7×7, 9×9, 11×11)作为MN数据集窗口获取的实验结果。通过获得的总体准确性(OA)，平均准确性(average accuracy, AA)和Kappa系数来分析窗口大小的影响。表 2显示了分类准确度随窗口大小增加的趋势，并且在9×9时获得了最佳性能。这意味着不合适的窗口大小无法合理地描述免疫复合物的分布特征。因此，在以下实验中，选择9×9作为MN数据集的窗口大小。

因实验数据为三阶张量，令降低维度为$\{ {B_1}, {B_2}, {B_3}\} $。与Pan等人(2017a)实验一致，空间维的降维设置为{1, 1}。通过固定窗口大小，$\lambda $，$\beta $和$\gamma $的值来调整光谱尺寸。图 4显示了不同降维方法(即NWFE(Kuo和Landgrebe，2004)，EFDC (Gao等，2012)，LFDA (Li等，2012)，BCGDA(Ly等，2014)，CDME(Lyu等，2017)，TSLGDA(Pan等，2017a))将MN数据集降维到不同维度后，采用SVM分类器对低维特征进行分类的结果。其中，Original对应未降维数据使用SVM分类获取的分类准确性。图 4表明大多数降维方法的OA随光谱维数的增加而增加，然后变得稳定。此外，与基于向量表示的降维方法(即NWFE，EFDC，LFDA，BCGDA，CDME)相比，基于张量表示的降维方法(即TSLGDA，T-LapSLRG)可以利用更多的空间结构信息并获得更好的结果。在大多数情况下，除了光谱维数为1的情况外，T-LapSLRG优于其他降维方法。下列实验使用的空间维度为{1, 1, 10}。

为了验证所提出的T-LapSLRG的有效性，采用SVM分类器(Moughal，2013)对获取的低维特征进行分类。将其分类性能与几种降维方法获得的低维特征通过SVM获取的分类性能进行比较，以证明T-LapSLRG的有效性。4个质量指标(即各个类别的准确性、总体准确性(OA)、平均准确性(AA)和Kappa系数(Kappa))用于衡量分类性能。表 3列出了通过所有降维方法获得的10维特征的OA和AA，这表明T-LapSLRG在这两个类别上均提供了最佳结果。详细而言，T-LapSLRG的性能优于其他降维方法，OA高出1.40 % ~34.75 %，AA高出1.46 % ~36.89 %，Kappa高出0.031~0.73。此外，结果证实基于张量表示的降维方法性能优于基于向量表示的降维方法。

表 4—表 6中列出了更详细的分类准确性评估。表 4列出了所有降维方法获得的HBV-MN患者每幅图像的OA。表 5列出了通过不同的降维方法获得的来自PMN患者的每幅图像的OA。从表 4和表 5可以看出，T-LapSLRG获得的性能优于其他方法。对于图像17002-3和17136-3，T-LapSLRG的OA小于85 %，但并不影响这两名患者的鉴定。表 6列出了通过不同方法获得的多幅图像的每个患者的OA。每个患者的OA指的是患者所有图像中所有像素的分类精度。例如，患者17002包括3幅HSI图像，分类结果则基于这3幅测试图像中的所有像素获得。表 6表明对于大多数患者而言，T-LapSLRG的识别准确性优于其他降维方法。此外，通过T-LapSLRG获得的所有患者的分类准确性均达到90 % 以上。在临床诊断中，可以在像素水平精度达到85 % 或更高时确定疾病的类型。因此得出，对于MN实验数据集，T-LapSLRG能够有效识别HBV-MN和PMN。

在实际医学应用中，可用训练样本的数量通常是有限的，因此研究算法对训练样本量的敏感度至关重要。图 5展示了通过不同降维方法(即NWFE，EFDC，LFDA，BCGDA，CDME，TSLGDA，T-LapSLRG)获取的低维特征通过SVM分类器获得的OA，训练样本的百分比从5 % 至50 % 不等。实验结果表明，T-LapSLRG始终比其他降维方法更有助于提升分类精度，尤其是在训练样本规模非常小(例如5 %)的情况下。综上，所有实验结果都证明了T-LapSLRG在实际MN临床鉴定中具有不可忽视的应用潜力。

本文提出了一种基于张量的拉普拉斯稀疏低秩图(T-LapSLRG)嵌入方法，并结合SVM分类器对医学高光谱图像进行判别分析。医学高光谱图像立方体被视为三阶张量，而基于张量的图则是通过从立方体中提取局部块(子张量)构建的。通过在目标函数上施加稀疏、低秩约束，T-LapSLRG可以同时挖掘局部和全局结构。通过引入流形正则项，T-LapSLRG可以调整图结构并进一步增强所获得的低维特征的判别能力。膜性肾病显微高光谱图像的实验结果表明，除了低维空间维度为1时，T-LapSLRG得到的分类精度都优于其他对比降维方法。当低维空间维度为10时，T-LapSLRG降维后获取的特征用SVM分类器对于PMN的分类精度为98.13 %，HBV-MN分类精度为96.60 %，OA为97.14 %，AA为97.05 %，Kappa为0.942。当训练样本比例从5 % 至50 % 变化时，本文算法得到的分类性能始终优于其他降维方法。即使对于较小的训练样本数量，T-LapSLRG方法也是有效的。实验结果证明，本文工作为膜性肾病诊断提供了一种新技术，具有临床应用潜力。然而，T-LapSLRG求解过程中使用的交替迭代算法时间复杂度较高，未来工作将集中于研究针对张量运算的加速算法。此外，本文算法在样本获取过程中依赖医生标注，这个过程会耗费大量人力成本。未来研究工作也将涉及部分医学样本辅助标记。