无线传感器网络（wireless sensor network，WSN）广泛应用于目标跟踪、环境监测等领域^[1]，其应用前提是能够提供精确的定位信息. 由于资源有限，WSN中只有少量信标节点能够通过配置GPS设备实现自身定位. 在这种情况下，可以利用信标节点的先验位置坐标及节点与信标节点之间的物理测度来实现对未知节点的定位. 目前有2类无线传感器网络定位技术^[2]，基于距离的定位通过使用不同的测距方案获得距离或角度信息，如信号接收强度（received signal strength，RSS）、到达时间（time of arrival，TOA）；基于距离无关的定位利用未知节点和信标节点之间的连通信息实现粗粒度的节点定位.

本研究主要针对基于距离的定位技术. 在实际情况下，该技术面临以下2个方面的挑战：1）由于环境和节点能量限制等因素的影响，部分节点对之间的距离信息缺失；2）在实际测距过程中，测距精度会受到噪声的影响. 测量过程中的噪声通常由高斯噪声、野值噪声和脉冲噪声组成. 其中，高斯噪声由硬件的局限性而产生，服从高斯分布；野值噪声由硬件故障、多径效应或恶意攻击等因素所引起，服从拉普拉斯分布，通常会导致严重的定位误差；脉冲噪声是环境的不确定性以及传感器节点信号收发模块故障所导致的欧氏距离矩阵(Euclidean distance matrix, EDM)中的连续误差，以行或列中部分连续误差的形式出现，服从拉普拉斯分布，其中连续误差的数量称为脉冲噪声的宽度. 通过以上分析可以看出，在实际情况下，由节点对之间距离信息构建出的距离矩阵（称为观测矩阵）是存在数据缺失的，并且是受到复合噪声污染的，无法直接用于节点定位. 因此，许多基于矩阵补全的算法^[3-4]被相继提出，其中部分算法^[4-5]考虑了噪声的影响，但是现有算法均忽略了脉冲噪声，致使定位精度还有待进一步提高.

针对上述问题，提出数据缺失和复合噪声条件下基于距离测度的无线传感器网络定位算法. 该算法分为2个阶段：1）以EDM的自然低秩性为基础，采用正则化技术，将EDM的重构问题转化为复合噪声条件下的矩阵补全问题，采用L_1,2范数显式平滑脉冲噪声，建立正则化矩阵补全模型. 为了有效求解该问题，将向量空间线性Bregman迭代算法推广到多维空间；2）在此基础上采用多维标度法（multi-dimensional scaling，MDS），设计鲁棒较高的基于Bregman散度的无线传感器网络定位算法（WSN localization via Bregman divergence based on MDS，LBD_MDS）.

WSN中的定位问题作为许多应用的关键部分，一直以来受到研究者们的格外关注. 定位精度是WSN定位问题面临的主要挑战之一，研究者们先后提出基于MDS^[6-7]、极大似然（maximum likelihood，ML）^[8]、指纹^[9-10]和半定规划（semi-definite programming，SDP）^[11]等的定位方法.

基于MDS的定位方法是在利用EDM生成节点相对坐标的基础上，对齐信标节点坐标，进而将节点相对坐标映射为绝对坐标^[6]，该方法对EDM的精度要求较高. Bhaskar^[8]将节点定位问题描述为概率问题，提出基于约束极大似然估计的算法，用于重建d维欧氏空间中的节点位置. Wang等^[9]研究RSS的时间相关性与定位精度之间的关系，并从理论上证明可以利用RSS的时间相关性来提高定位的准确性. Singh等^[12]针对分布式、各向同性的WSN，仅使用单个锚节点作为参考节点，提出利用虚拟锚节点投影来实现节点定位的概念，在一定程度上解决了定位过程中的视距遮挡问题. 上述方法均以相对精确的节点间距离测度为前提实现节点定位，对于存在数据缺失和复合噪声污染的EDM，这些方法不适用.

针对存在数据缺失和复合噪声污染的EDM，Arora等^[13]引入蝶形优化算法解决高斯噪声干扰下的无线传感器网络定位问题. Fang等^[7]利用自适应卡尔曼滤波消除测量噪声的影响，提出基于MDS和自适应卡尔曼滤波的定位算法，以较高的定位精度和较低的时间复杂度实现节点定位. 另外，他们针对噪声环境中加权k近邻算法不能较好地进行节点位置估计的情况^[10]，利用自适应卡尔曼滤波和模因算法，提出噪声条件下WSN指纹定位的最优加权k近邻算法. Sahota等^[14]考虑多路径效应的影响，给出网络中各种多媒体和多径通信场景的路径丢失和衰落模型，根据传输距离和节点的位置坐标对接收信号强度进行建模，在此基础上，基于极大似然优化，利用导出的统计模型实现节点定位. 上述研究均力求在定位过程中减小噪声对定位精度的影响，而未涉及到距离测度中的噪声检测和分离. Feng等^[3]将矩阵补全理论引入无线传感器定位，将定位问题转化为低秩矩阵恢复问题，但其仅简单地将高斯噪声作为测量噪声，未考虑复合噪声的情况，进而导致定位精度不高. Guo等^[11]考虑测量噪声和异常值对EDM的影响，利用半定嵌入定理，提出基于SDP的低秩矩阵补全算法，并推导加权半定松弛定位方法，提高节点定位的精度. 但是该算法复杂度较高，不适用于处理大规模WSN定位问题. Xiao等^{[4, 15]}将高斯噪声和野值点噪声作为复合噪声同时考虑，基于范数正则化技术和交替方向乘子法，实现EDM的重构，虽然定位精度有一定程度的提高，但忽略了脉冲噪声的影响，致使定位结果仍不理想. 针对上述情况，本研究采用正则化技术构建矩阵补全模型，并基于扩展的线性Bregman迭代方法，设计实现鲁棒较高的WSN节点定位算法，以提高复合噪声干扰下的定位精度.

线性Bregman迭代作为一种优化算法被广泛应用于压缩感知^[16]、图像去噪^[17]、目标检测^[18]和量化聚类^[19]等领域，已经成为求解范数优化问题最有效的方法之一.

定义1 　Bregman散度^[20]. 令函数 $J\left( { x} \right):{{\bf R}^n} \to {\bf R}$为连续的凸函数，任给函数上的两点 ${ u}$、 ${ v}$（ $\forall { u},{ v} \in { x}$），函数 $J$在 ${ u}$、 ${ v}$两点之间的Bregman散度表达式为

(1) $ D_J^{ p}\left( {{ u},{ v}} \right) = J\left( { u} \right) - J\left( { v} \right) - \left\langle {{ p},{ u} - { v}} \right\rangle . $

式中： $\left\langle { \cdot , \cdot } \right\rangle $为内积， ${ p} \in \partial J\left( { v} \right)$为函数 $J$在 ${ v}$点的次梯度， $\partial J\left( { v} \right)$为函数在 ${ v}$点的次微分.

定义2 　多元函数Bregman散度.令多元函数 $J\left( {{{ x}_1},{{ x}_2}, \cdots ,{{ x}_\ell }} \right):{{\bf R}^n} \to {\bf R}$为连续的凸函数， $\forall {{ u}_i},{{ v}_i} \in $ $ {{ x}_i},i = 1,2, \cdots ,\ell $，多元函数 $J$ 在 $( {{{ u}_1},}{{{ u}_2}, \cdots ,{{ u}_\ell }} )$、 $( {{ v}_1}, $ ${{ v}_2}, \cdots ,{{ v}_\ell })$两点之间的Bregman散度表达式为

(2) $ \begin{split} &\!\! D_J^{ P}\left[ {\left( {{{ u}_1},{{ u}_2}, \cdots ,{{ u}_\ell }} \right),\left( {{{ v}_1},{{ v}_2}, \cdots ,{{ v}_\ell }} \right)} \right] = \\ &\!\! J\left( {{{ u}_1},{{ u}_2}, \cdots ,{{ u}_\ell }} \right) - J\left( {{{ v}_1},{{ v}_2}, \cdots ,{{ v}_\ell }} \right) - \displaystyle\sum\limits_{i = 1}^\ell {\left\langle {{{ p}_{{{ v}_i}}},{{ u}_i} - {{ v}_i}} \right\rangle .} \end{split} $

式中： ${ P} = \left( {{{ p}_{{{ v}_1}}},{{ p}_{{{ v}_2}}}, \cdots, {{ p}_{{{ v}_\ell }}}} \right) \in \partial J\left( {{{ v}_1},{{ v}_2}, \cdots ,{{ v}_\ell }} \right)$为函数 $J$在 $\left( {{{ v}_1},{{ v}_2}, \cdots ,{{ v}_\ell }} \right)$点的次梯度.

以下为3个函数的例子，并给出相应的Bregman散度（ $D_J^{ P}\left[ {\left( {{{ u}_1},{{ u}_2}, \cdots ,{{ u}_\ell }} \right),\;\left( {{{ v}_1},{{ v}_2}, \cdots ,{{ v}_\ell }} \right)} \right]$，简记为 $D_J^{ P}\left( {{ U},{ V}} \right)$）.

(3) $ \left. \begin{array}{l} J\left( {{{ x}_1},{{ x}_2}, \cdots ,{{ x}_\ell }} \right) = \displaystyle\sum\limits_{i = 1}^\ell {{{\left\| {{{ x}_{{i}}}} \right\|}^2}}, \\ D_J^{ P}\left( {{ U},{ V}} \right) = \displaystyle\sum\limits_{i = 1}^\ell {{{\left\| {{{ u}_i}} \right\|}^2}} - \displaystyle\sum\limits_{i = 1}^\ell {{{\left\| {{{ v}_i}} \right\|}^2}} - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; \displaystyle\sum\limits_{i = 1}^\ell {\left\langle {2{{ v}_i},{{ u}_i} - {{ v}_i}} \right\rangle } = \displaystyle\sum\limits_{i = 1}^\ell {{{\left\| {{{ u}_i} - {{ v}_i}} \right\|}^2}} . \\ \end{array} \right\} $

当 $\ell = 1$时， $D_J^P\left( {{ U},{ V}} \right)$为欧氏距离的平方.

(4) $ \left. \begin{array}{l} J\left( {{{ x}_1},{{ x}_2}, \cdots ,{{ x}_\ell }} \right) = \displaystyle\sum\limits_{i = 1}^\ell {{{ x}_{{i}}}^{\rm T}{ A}} {{ x}_i},\\ D_J^{ P}\left( {{ U},{ V}} \right) = \displaystyle\sum\limits_{i = 1}^\ell {{{\left( {{{ u}_i} - {{ v}_i}} \right)}^{\rm T}}{ A}} \left( {{{ u}_i} - {{ v}_i}} \right). \end{array} \right\} $

当 $\ell = 1$时， $D_J^{ P}\left( {{ U},{ V}} \right)$为马氏距离.

(5) $ \left. \begin{array}{l} J\left( {{{ x}_1},{{ x}_2}, \cdots ,{{ x}_\ell }} \right) = \displaystyle\sum\limits_{i = 1}^\ell {\displaystyle\sum\limits_{j = 1}^n {\left( {{{x}_{{{i}}j}}\ln \;{{x}_{ij}}} \right)} } ,\\ D_J^{ P}\left( {{ U},{ V}} \right) = \displaystyle\sum\limits_{i = 1}^\ell {\displaystyle\sum\limits_{j = 1}^n {\left( {{{u}_{{{i}}j}}\ln\; \displaystyle\frac{{{{u}_{ij}}}}{{{{v}_{ij}}}}} \right)} } . \end{array} \right\} $

当 $\ell = 1$时， $D_J^{ P}\left( {{ U},{ V}} \right)$为KL散度.

对于部署于 $d$维空间的某区域 $S$内（ $S \in {{\bf R}^d}$）的无线传感器网络，假设 $n$个无线传感器节点随机部署于 $S$内，设 ${ X} = \left[ {{{ x}_1},{{ x}_2}, \cdots ,{{ x}_n}} \right]$（ ${ X} \in {{\bf R}^{d \times n}}$，为 $d$维空间中 $n$个节点的坐标矩阵），可以求得欧氏距离矩阵 ${ R} \in {{\bf R}^{n \times n}}$（其中 ${R_{ij}} = \; {\left\| {{{ x}_i} - {{ x}_j}} \right\|^2}$ $,\;i,j = 1,2, \cdots ,n$）. 构建节点间EDM的观测矩阵 ${ M} \in {{\bf R}^{n \times n}}$，如上所述，矩阵 ${ M}$是不完整的、含有噪声的. 将无线传感器网络定位分为2个阶段：EDM恢复、节点坐标映射. 由于矩阵 ${ M}$的不完整性和含噪特性，对 ${ M}$简单地采用基于MDS的定位方法不能对无线传感器网络实现精确定位，须基于矩阵补全技术利用观测矩阵 ${ M}$得到精确的欧氏距离矩阵 ${ R}$的估计.

Fu等^[21]已经给出 ${\rm{rank}}\left( { R} \right) \leqslant d + 2$的证明（ ${\rm{rank}}\left( { R} \right)$为矩阵 ${ R}$的秩），因此在 $n \gg d$的情况下，矩阵 ${ R}$是低秩的. EDM恢复问题可以建模为如下矩阵补全模型：

(6) $\left. \begin{array}{l} \mathop {\min \; }\limits_{{ R},{ G},{ O},{ C} \in {\mathbb{R}^{n \times n}}} \;\;\;\;\;\; {\left\| { R} \right\|_*} + \varphi \left\| { G} \right\|_{\rm F}^2 + \mu {\left\| { O} \right\|_1} + \lambda {\left\| { C} \right\|_{1,2}}. \\ {\rm{s.t.}}\;\;\;\;\; {P_{\varOmega} }\left( { M} \right) = {P_{\varOmega} }\left( {{ R} + { G} + { O} + { C}} \right). \\ \end{array}\right\} $

式中： $ { G}{\text{、}}{ O}{\text{、}}{ C}$分别为高斯噪声矩阵、野值噪声矩阵和脉冲噪声矩阵，其中脉冲噪声包括行脉冲噪声和列脉冲噪声； ${\left\|{ R}\right\|}_*{\text{、}}{\left\|{ G} \right\|}_{\rm F}{\text{、}}\left\|{ O} \right\|_1{\text{、}}\left\|{ C} \right\|_{1,2}$分别为矩阵核范数、Frobenius范数、L₁范数、L_1,2范数； $\varphi {\text{、}}\mu {\text{、}}\lambda $为平衡3种噪声的可调参数； $\varOmega \subseteq \left[ {1:n} \right] \times $ $ \left[ {1:n} \right]$为观测元素的索引集合，若 $\left( {i,j} \right) \in \varOmega $，则 ${\left[ {{P_{\varOmega} }\left( { M} \right)} \right]_{ij}}{\rm{ = }}{ M}_{ij}$；否则 ${\left[ {{P_{\varOmega} }\left( { M} \right)} \right]_{ij}}{\rm{ = }}{0}$. 对式（6）进行求解，可得到欧氏距离矩阵 ${ R}$的估计 ${\hat { R}}$（ ${\hat { R}} \in {{\bf R}^{n \times n}}$）.

引入Bregman散度求解式（6）中的矩阵补全模型. 式（6）对应的增广拉格朗日函数表达式为

(7) $ \begin{array}{l} {L_\rho }\left( {{ R},{ G},{ O},{ C},{ Y}} \right)\; = \; {\left\| { R} \right\|_*} + \varphi \left\| { G} \right\|_{\rm F}^2 + \mu {\left\| { O} \right\|_1} + \lambda {\left\| { C} \right\|_{1,2}} + \\ \;\;\;\; \left\langle {{ Y},{ R} + { G} + { O} + { C} - { M}} \right\rangle + {{\rho \left\| {{ M} - \left( {{ R} + { G} + { O} + { C}} \right)} \right\|_{\rm F}^2} / 2}. \\ \end{array} $

式中： ${ Y} \in {{\bf R}^{n \times n}}$为拉格朗日乘子； $\rho > 0$为可调参数，大小与高斯噪声项呈负相关，若将 $\rho $设置为较大的值，可达到隐式平滑高斯噪声的目的. 因此可以将式（6）简化为如下模型：

(8) $\left. \begin{array}{l} \mathop {\min \; }\limits_{{ R},{ O},{ C} \in {{\bf R}^{n \times n}}} \;\;\; {\left\| { R} \right\|_*} + \mu {\left\| { O} \right\|_1} + \lambda {\left\| { C} \right\|_{1,2}}. \\ \;\;\;\;\; {\rm{s.t.}}\;\;\;\;\; {P_{\varOmega} }\left( { M} \right) = {P_{\varOmega} }\left( {{ R} + { O} + { C}} \right). \end{array}\right\} $

为了有效求解式（8），将其松弛为如下无约束优化问题：

(9) $\begin{split} \mathop {\min \; }\limits_{{ R},{ O},{ C} \in {{\bf R}^{n \times n}}} &\;\; \tau \;\; \left( {{{\left\| { R} \right\|}_*} + \mu {{\left\| { O} \right\|}_1} + \lambda {{\left\| { C} \right\|}_{1,2}}} \right){\rm{ + }} \\ & {{\left\| {{P_{\varOmega} }\left( {{ M} - \left( {{ R} + { O} + { C}} \right)} \right)} \right\|_{\rm F}^2} / 2} . \end{split} $

式中： $ \tau > 0$.

为了描述方便，令

(10) $\left. \begin{array}{l} J\left( {{ R},{ O},{ C}} \right) = \tau \;\; \left( {{{\left\| { R} \right\|}_*} + \mu {{\left\| { O} \right\|}_1} + \lambda {{\left\| { C} \right\|}_{1,2}}} \right), \\ H\left( {{ R},{ O},{ C}} \right) = {{\left\| {{P_{\varOmega} }\left( {{ M} - \left( {{ R} + { O} + { C}} \right)} \right)} \right\|_{\rm F}^2} / 2}. \end{array} \right\} $

则式（9）可改写为

(11) $ \mathop {\min \; }\limits_{{{ R},{ O},{ C}} \in {{\bf R}^{n \times n}}} J\left( {{ R},{ O},{ C}} \right) + H\left( {{ R},{ O},{ C}} \right). $

引入多元函数Bregman散度求解式（11）. 根据定义2，函数 $J$在 $\left( {{ R},{ O},{ C}} \right)$、 $\left( {{{ R}^k},{{ O}^k},{{ C}^k}} \right)$两点间的Bregman散度表达式为

(12) $ \begin{split} & D_{ J}^{{{ P}^k}}\left[ {\left( {{ R},{ O},{ C}} \right),\left( {{{ R}^k},{{ O}^k},{{ C}^k}} \right)} \right] = \\ & \;\;\;\quad J\left( {{ R},{ O},{ C}} \right) - J\left( {{{ R}^k},{{ O}^k},{{ C}^k}} \right) - A. \end{split} $

式中： $A = \left\langle {\tau { P}_{ R}^k,{ R} - {{ R}^k}} \right\rangle {\rm{ + }}\left\langle {\tau \mu { P}_{ O}^k,{ O} - {{ O}^k}} \right\rangle + $ $\left\langle {\tau \lambda { P}_{ C}^k},\right.{ C} -$ $\left. {{ C}^k} \right\rangle , $ ${{ P}^k} = \; \left( {{ P}_{ R}^k,{ P}_{ O}^k,{ P}_{ C}^k} \right) \in \partial J\left( {{{ R}^k},{{ O}^k},{{ C}^k}} \right)$为函数 $J$在 $ ( {{ R}^k},$ ${{ O}^k},{{ C}^k} )$点的次梯度.

令 $E = \left( {{ R},{ O},{ C}} \right)$， ${E^k} =\left( {{{ R}^k},{{ O}^k},{{ C}^k}} \right)$，式（11）可按照如下方式进行迭代求解：

(13) $\left. \begin{aligned} & \left( {{{ R}^{k + 1}},{{ O}^{k + 1}},{{ C}^{k + 1}}} \right){\rm{ = }} \mathop {\arg \min }\limits_{{ R},{ O},{ C} \in {{\bf R}^{n \times n}}} D_J^{{{ P}^k}}\left( {E,{E^k}} \right) + H\left( E \right) . \\ {\rm{s.t.}}&\left\{ \begin{aligned} 0 \in \partial \left( {D_J^{{{ P}^k}}\left( {E,{E^k}} \right) + H\left( E \right)} \right)\left| {_{{{ R}^{k + 1}}}} \right. ; \\ 0 \in \partial \left( {D_J^{{{ P}^k}}\left( {E,{E^k}} \right) + H\left( E \right)} \right)\left| {_{{{ O}^{k + 1}}}} \right. ; \\ 0 \in \partial \left( {D_J^{{{ P}^k}}\left( {E,{E^k}} \right) + H\left( E \right)} \right)\left| {_{{{ C}^{k + 1}}}} \right. . \end{aligned} \right.\end{aligned} \right\} $

受分裂Bregman迭代^[22]思想启发，将其应用于矩阵空间，同时引入交替最小化方法求解 $\left( {{{ R}^{k + 1}},{{ O}^{k + 1}},{{ C}^{k + 1}}} \right)$，则式（11）可以进一步采用算法1中所示的分裂Bregman迭代-交替最小化算法（split Bregman iteration- alternating minimization, SBI-AM）进行求解. 其中， $N$为最大迭代次数， $({{ R}^{\rm{opt}}},{{ O}^{\rm{opt}}},{{ C}^{\rm{opt}}}$)为算法的解.

算法1　SBI-AM算法

输入： ${{ O}^0} = {{ C}^0}{\rm{ = {\bf 0}}}$， $P_{ O}^0 = P_{ C}^0 = {\bf 0}$， ${P_{\varOmega} }\left( {{ M}} \right)$，N

输出： ${{ R}^{\rm{opt}}},{{ O}^{\rm{opt}}},{{ C}^{\rm{opt}}}$

1. for k=0 to N

2. ${{ R}^{k + 1}} = \mathop {\arg \min }\limits_{{ R} \in {{\bf R}^{n \times n}}}\; \tau {\left\| { R} \right\|_*} - \tau \left\langle {{ P}_{ R}^k,{ R}} \right\rangle + {{\left\| {{P_{\varOmega} }}\right.}}({ M} - { R} -$ $\left. \; { O} - { C} ) \right\|_{\rm F}^2 /2 $

3. $ {{ O}^{k\! +\! 1}} \!\!=\!\! \mathop {\arg \min }\limits_{{ O} \in {{\bf R}^{n \times n}}}\;\! \tau \mu {\left\| { O} \right\|_1}\!\! -\!\! \tau \mu \left\langle {{ P}_{ O}^k,{ O}} \right\rangle \! +\! \left\| P_{\varOmega}\left( { M}\!-\!{ { R}^{k \!+\! 1}\! -\! }{ O}\! -\!{ C} \right) \right\|_{\rm F}^2 \!\bigg/ \!2$

$4.\; {{ C}^{k + 1}} = \mathop {\arg \min }\limits_{{ C} \in {{\bf R}^{n \times n}}}\; \tau \lambda \;{\left\| { C} \right\|_{1,2}} - \tau \lambda\;\left\langle {{ P}_{ C}^k,{ C}} \right\rangle \; + \; \left\| P_{\varOmega}\;\left( { M}-{ {{ R}^{k + 1}} - } \quad\quad\right.\right. $ $\left. \left.{ O}^{k+1} -{ C} \right) \right\|_{\rm F}^2 \Big/ 2 $

5. ${ P}_{ R}^{k + 1} = { P}_{ R}^k + {{{P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^{k + 1}} - {{ C}^{k + 1}}} \right)} / \tau }$

6. ${ P}_{ O}^{k + 1} = { P}_{ O}^k + {{{P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^{k + 1}} - {{ C}^{k + 1}}} \right)} / {(\tau \mu) }}$

7. ${ P}_{ C}^{k + 1} = { P}_{ C}^k + {{{P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^{k + 1}} - {{ C}^{k + 1}}} \right)} / {(\tau \lambda) }}$

8. end for

9. return ${{ R}^{\rm{opt}}} \leftarrow {{ R}^{N + 1}},{{ O}^{\rm{opt}}} \leftarrow {{ O}^{N + 1}},{{ C}^{\rm{opt}}} \leftarrow {{ C}^{N + 1}}$

由SBI-AM算法不难看出，由于函数 ${\left\| { R} \right\|_*}$、 ${\left\| { O} \right\|_1}$、 ${\left\| { C} \right\|_{1,2}}$是不可微的，算法中的第2~4步不能直接求解相应变量. 为了进一步求解变量，给出如下基本定义和定理：

定义3 　近邻算子^[23-24].令 $g\left( { X} \right)$为 ${{\bf R}^{m \times n}}$上的实值凸函数， $\tau > 0$， $\forall { Z} \in {{\bf R}^{m \times n}}$， $g\left( { X} \right)$的近邻算子表达式为

(14) $ {{\rm{prox}}_{\tau g\left( { X} \right)}}\left( { Z} \right)\; = \; \mathop {\arg \min }\limits_{{ X} \in {{\bf R}^{m \times n}}} \;\left( {\tau g\left( { X} \right) + {{\left\| {{ X} - { Z}} \right\|_{\rm F}^2} / 2}} \right). $

定理1^[25] 　设 $ F_1\text{、}F_2$为 ${{\bf R}^{m \times n}}$上2个下半连续的凸函数，且 ${F_2}$在 ${{\bf R}^{m \times n}}$中可微并具有β-Lipschitz连续梯度，则对于凸优化问题 $\mathop {\min }\limits_{{ X} \in {{\bf R}^{m \times n}}} {F_1}\left( { X} \right) + {F_2}\left( { X} \right)$，若函数 ${F_1} + {F_2}$是强制的且严格凸的，则该问题存在唯一解，且对任意初始值 ${{ X}_0}$及 $0 < \delta < {2 / \beta }$，采用如下方法生成的迭代序列 ${{ X}_{k + 1}}$能够收敛到该问题的唯一解：

(15) $ {{ X}_{k + 1}} = \mathop {\arg \min }\limits_{{ X} \in {{\bf R}^{m \times n}}}\; \left( {\delta {F_1}\left( { X} \right) + \;\; {{\left\| {{ X}\; - \; \left( {{{ X}_k}\; - \; \delta \nabla {F_2}\left( {{{ X}_k}} \right)} \right)} \right\|_{\rm F}^2} / 2}} \right). $

定理2^[26] 　对于任意的 $\kappa > 0$， ${ Z} \in {{\bf R}^{m \times n}}$，矩阵 ${ X}$的核范数的近邻算子表达式为

(16) $ {{\rm{prox}}_{\kappa {{\left\| { X} \right\|}_*}}}\left( { Z} \right) = {D_\kappa }\left( { Z} \right). $

式中： ${D_\kappa }\left( \cdot \right)$为奇异值阈值算子^[26].

定理3^[24] 　对任意 $\kappa > 0$， ${ Z} \in {{\bf R}^{m \times n}}$，矩阵 ${ X}$的 ${L_1}$范数的近邻算子表达式为

(17) $ {{\rm{prox}}_{\kappa {{\left\| { X} \right\|}_1}}}\left( { Z} \right) = {S_\kappa }\left( { Z} \right). $

式中： ${S_\kappa }\left( \cdot \right)$为收缩算子^[26].

定理4^[27] 　对任意 $\kappa > 0$， ${ X},{ Z} \in {{\bf R}^{m \times n}}$，矩阵 ${ X}$的 ${L_{1,2}}$范数的近邻算子表达式为

(18) $ {{\rm{prox}}_{\kappa {{\left\| { X} \right\|}_{1,2}}}}\left( { Z} \right) = {J_\kappa }\left( { Z} \right). $

式中： ${J_\kappa }\left( {{ Z}} \right){_{\left( i \right)}} \!=\! \max\, \left\{ {1 \!-\! {{\kappa } \Big/ {{{\left\| {{{ Z}_{\left( i \right)}}} \right\|}_2},0}}} \right\}{{ Z}_{\left( i \right)}},$ $i = 1,2, \cdots ,m$.

采用上述定义和定理，对SBI-AM算法中的变量 $\left( {{{ R}^{k + 1}},{{ O}^{k + 1}},{{ C}^{k + 1}}} \right)$进行求解，步骤如下.

1）更新R.

令 $\varGamma _{ R}^k = \delta {P_{\varOmega} }\left( {{ M} - {{ R}^k} - {{ O}^k} - {{ C}^k}} \right)$，由定义3和定理1， ${{ R}^{k + 1}}$、 ${ P}_{ R}^{k + 1}$被修正为

(19) $\left.\begin{split} &\!\!\!\!\!\!\!\!{{ R}^{k + 1}}\! =\! \mathop {\arg \min }\limits_{ R}\!\left[ \tau \delta {\left\| { R} \right\|_*} \!+\! {{\left\| {{ R} \!-\! {{ R}^k}\! -\! \tau \delta { P}_{ R}^k - \varGamma _{ R}^k} \right\|_{\rm F}^2} \Big/ 2}\right],\\ &\!\!\!\!\!\!\!\!{ P}_{ R}^{k + 1} = { P}_{ R}^k - {{\left( {{{ R}^{k + 1}} - {{ R}^k} - \varGamma _{ R}^k} \right)} \Big/ {(\tau \delta) }}. \end{split}\right\}$

令

(20) $ {{ B}^k} = {{ R}^k} + \tau \delta { P}_{ R}^k + \varGamma _{ R}^k\, , $

有 $ {{ R}^{k + 1}} = \mathop {\arg \min }\limits_{ R}\;\left[ \tau \delta {\left\| { R} \right\|_*} + {{\left\| {{ R} - {{ B}^k}} \right\|_{\rm F}^2} \Big/ 2}\right].$

根据定理2可以得到

(21) $ {{ R}^{k + 1}} = {D_{\tau \delta }}\left( {{{ B}^k}} \right). $

对于 $ { B}$，有

(22) $ {{ B}^k} - {{ B}^{k - 1}} = \delta {P_{\varOmega} }\left( {{ M} - {{ R}^k} - {{ O}^k} - {{ C}^k}} \right). $

2）更新O. 与步骤1）相似，对于矩阵 ${ O}$，

(23) $ {{ O}^{k + 1}} = \mathop {\; \arg \min }\limits_{ O} \;\left[\tau \mu \delta {\left\| { O} \right\|_1} + \; {{\left\| {{ O} - {{ O}^k} - \tau \mu \delta { P}_{ O}^k - \varGamma _{ O}^k} \right\|_{\rm F}^2}\Big/ 2}\right]. $

式中： $ \varGamma _{ O}^k = \delta {P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^k} - {{ C}^k}} \right)$.

令 $ {{ I}^k} = {{ O}^k} + \tau \mu \delta { P}_{ O}^k + \varGamma _{ O}^k$，

有

(24) $\left. \begin{aligned} &{{ O}^{k + 1}} = \mathop {\arg \min }\limits_{ O}\;\left[ \tau \mu \delta {\left\| { O} \right\|_1} + {{\left\| {{ O} - {{ I}^k}} \right\|_{\rm F}^2} \Big/ 2}\right],\\ &{{ I}^k} - {{ I}^{k - 1}} = \delta {P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^k} - {{ C}^k}} \right). \end{aligned} \right\} $

由定理3，可以得到式（24）的解析解：

(25) $ {{ O}^{k + 1}} = {S_{\tau \mu \delta }}\left( {{{ I}^k}} \right). $

3）更新C. 同理，对于脉冲噪声矩阵 ${ C}$，

(26) $ {{ C}^{k + 1}} = \mathop {\arg \min }\limits_{ C}\;\left[ \tau \lambda \delta {\left\| { C} \right\|_{1,2}} + \;\; {{\left\| {{ C} - {{ C}^k} - \tau \lambda \delta { P}_{ C}^k - \varGamma _{ C}^k} \right\|_{\rm F}^2} \Big/2}\right]. $

式中： $ \varGamma _{ C}^k = \delta {P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^{k + 1}} - {{ C}^k}} \right)$.

令 ${{ U}^k} = {{ C}^k} + \tau \lambda \delta { P}_{ C}^k + \varGamma _{ C}^k$，可以得到

(27) $ \left. \begin{aligned} &{{ C}^{k + 1}} = \mathop {\arg \min }\limits_{ C}\;\left[ \tau \lambda \delta {\left\| { C} \right\|_{1,2}} + {{\left\| {{ C} - {{ U}^k}} \right\|_{\rm F}^2} \Big/ 2}\right],\\ &{{ U}^k} = {{ U}^{k - 1}} + \delta {P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^{k + 1}} - {{ C}^k}} \right). \end{aligned}\right\}$

由定理4，可以得到式（27）的解析解

(28) $ {{ C}^{k + 1}} = {J_{\tau \lambda \delta }}\left( {{{ U}^k}} \right). $

将上述步骤进行归纳，求解式（6）的整个优化过程可以整理为如算法2所示的基于Bregman散度的矩阵补全算法（matrix completion via Bregman divergence，BDMC）. 由算法2可以看出，BDMC算法是鲁棒高效的，且具有一定的可扩放性. 主要表现为：一方面， ${{ R}^{k + 1}}$的低秩性和 $\left( {{{ O}^{k + 1}},{{ C}^{k + 1}}} \right)$的稀疏性在迭代过程中能够得到保持；另一方面，对于复杂的矩阵奇异值分解，BDMC算法在每次迭代中仅涉及一次.

算法2　BDMC算法

输入： ${P_{\varOmega} }\left( { M} \right)$， $\tau ,\mu ,\lambda ,\delta $，N

输出： ${{ R}^{\rm{opt}}},{{ O}^{\rm{opt}}},{{ C}^{\rm{opt}}}$

1. Initialize ${{ R}^0} = {\bf 0},{{ O}^0} = {\bf 0},{{ C}^0} = {\bf 0},{{ B}^{ - 1}} = {\bf 0},$ ${{ I}^{ - 1}} = {\bf 0},{{ U}^{ - 1}} = {\bf 0}$

2. for k=0 to N

3. ${{ B}^k} = {{ B}^{k - 1}} + \delta {P_{\varOmega} }\left( {{ M} - {{ R}^k} - {{ O}^k} - {{ C}^k}} \right)$

4. ${{ R}^{k + 1}} = {D_{\tau \delta }}\left( {{{ B}^k}} \right)$

5. ${{ I}^k} = {{ I}^{k - 1}} + \delta {P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^k} - {{ C}^k}} \right)$

6. ${{ O}^{k + 1}} = {S_{\tau \mu \delta }}\left( {{{ I}^k}} \right)$

7. ${{ U}^k} = {{ U}^{k - 1}} + \delta {P_{\varOmega} }\left( {{ M} - {{ R}^{k + 1}} - {{ O}^{k + 1}} - {{ C}^k}} \right)$

8. ${{ C}^{k + 1}} = {J_{\tau \lambda \delta }}\left( {{{ U}^k}} \right)$

9. end for

10. return ${{ R}^{\rm{opt}}} \leftarrow {{ R}^{N + 1}},{{ O}^{\rm{opt}}} \leftarrow \; {{ O}^{N + 1}},{{ C}^{\rm{opt}}}\leftarrow $ $ \;\; {{ C}^{N + 1}} $

采用BDMC算法得到EDM的估计 ${\hat { R}}$，利用 ${\hat { R}}$中对应节点间的距离测度能够准确估算节点间的相对位置信息. 在此基础上，采用基于MDS的定位方法对WSN节点实现定位. 如算法3所示为基于Bregman散度的无线传感器网络定位算法的详细步骤. 其中， $k$为信标节点的个数， ${{ T}_i}$为第 $i$个未知节点的坐标， $i = k + 1,k + 2, \cdots ,n$.

算法3　LBD_MDS算法

输入： ${P_{\varOmega} }\left( { M} \right),$ $\tau ,\;\mu ,\;\lambda ,\;\delta, $ N，信标节点坐标 $ \{{T_1},{T_2}, \cdots ,$ ${T_k}| {k \geqslant 3} \}.$

输出：所有未知位置信息的节点坐标 $ \{ {T_i}| i = k + 1,$ $k + 2, \cdots ,n \}.$

　/^*重构EDM^*/

1. 采用BDMC算法，得到EDM的估计 ${{\hat R}} $

2. $\left[ {{ U},{ \varLambda} ,{ V}} \right] = {\rm{svd}}\left( { - {{{ \varTheta} { {\hat R}}{{ \varTheta} ^{\rm T}}} / 2}} \right)$

${ \varTheta} = { I} - {1 / {\rm{n}}}\left( {{\bf 1} \cdot {{\bf 1}^{\rm T}}} \right)$， ${ I}$为单位阵

3. 计算相对坐标矩阵

${ W} = { \varLambda} _d^{{1 / 2}}{ U}{\left( {:,1:d} \right)^{\rm T}}$

${ W} = \left[ {{{ W}_1},{{ W}_2}, \cdots ,{{ W}_n}} \right] \in {{\bf R}^{d \times n}}$

4. 计算变换矩阵

${ Q} = \displaystyle\frac{{\left[ {{{ T}_2} - {{ T}_1},{{ T}_3} - {{ T}_1}, \cdots, {{ T}_k} - {{ T}_1}} \right]}}{{\left[ {{{ W}_2} - {{ W}_1},{{ W}_3} - {{ W}_1}, \cdots, {{ W}_k} - {{ W}_1}} \right]}}$

5. 节点坐标映射

$\left\{ {{ T}\left| {{{ T}_i} - {{ T}_1} = { Q} \times \left( {{{ W}_i} - {{ W}_1}} \right),i = 1,2, \cdots ,n} \right.} \right\}$

6. return ${ T}$

为了评估LBD_MDS算法的性能，选取EDM恢复误差、定位误差、定位误差方差和定位误差累积分布等作为评估指标，与IALM^[28]、OptSpace^[29]、ScGrassMC^[30]进行比较. 假设WSN在边长为100 m的正方形区域内随机分布有100个节点，其中部分节点为信标节点，根据节点间的距离信息得到EDM. 向EDM中添加噪声，并以采样率sr对含噪EDM进行随机采样，得到观测矩阵 ${ M}$，以观测矩阵 ${ M}$作为上述算法的训练数据. 为了避免偶然性，重复进行20次试验，取平均值作为实验结果. 根据噪声环境的不同，设置如下4种情况：1）假设EDM未受到任何噪声污染，即观测矩阵 ${ M}$中的值为准确值. 2）假设EDM受到高斯噪声和野值噪声的污染. 高斯噪声服从均值为0，方差为100的高斯分布；野值噪声服从均值为0，方差为10 000的拉普拉斯分布. 3）假设EDM受到脉冲噪声的污染. 脉冲噪声宽度为30，服从均值为0，方差为10 000的拉普拉斯分布. 4）假设EDM受到高斯噪声、野值噪声和脉冲噪声的污染. 高斯噪声服从均值为0，方差为100的高斯分布；野值噪声服从方差为10 000，均值为0的拉普拉斯分布；脉冲噪声宽度为30，服从均值为0，方差为10 000的拉普拉斯分布.

令 ${ X} \in {{\bf R}^{2 \times n}}$， ${ R} \in {{\bf R}^{n \times n}}$（ $n = 100$为节点数量）分别为节点坐标矩阵和EDM. 选取以下4个指标评估LBD_MDS算法的性能.

1）EDM恢复误差：

(29) $ {\rm{REs}} = {{{{\left\| {{\hat { R}} - { R}} \right\|}_{\rm F}}} \Big/ {{{\left\| { R} \right\|}_{\rm F}}}}. $

式中： ${\hat { R}}$为通过矩阵补全算法得到的EDM的估计.

2）节点定位误差：

(30) $ {\rm{LEs}} = {{{{\left\| {{\hat { X}} - { X}} \right\|}_{\rm F}}} \Big/ n}. $

式中： ${\hat { X}}$为节点坐标矩阵 ${ X}$的估计.

3）节点定位误差方差：

(31) $ {\rm{LEV}} = {{\sum\nolimits_{i = 1}^n {{{\left( {\Delta {L_i} - {\rm{LEs}}} \right)}^2}} } / {{n}}}. $

式中： $\Delta {L_i} = \left[ {{{\left( {{{\hat x}_i} - {x_i}} \right)}^2} + {{\left( {{{\hat y}_i} - {y_i}} \right)}^2}}\right]^{1/2} $为第 $i$个节点的定位误差， $\left( {{x_i},{y_i}} \right)$、 $\left( {{{\hat x}_i},{{\hat y}_i}} \right)$（ $i = 1,2, \cdots ,n$）分别为第 $i$个节点坐标及其估计.

4）节点定位误差累积分布：

(32) $ {\rm{LE\_CDF}} = P_{ {\Delta {L_i} \leqslant \sigma } }. $

式中： $ \sigma$为常数.

实验1 　算法收敛速度比较. 在采样率sr为30%、50%的条件下，各算法的收敛情况如图1所示. 图中，Nc为算法的迭代次数，最大迭代次数设置为2 000，恢复误差阈值设置为10⁻⁸. 可以看出，相对于另外3种算法，ScGrassMC算法收敛速度最快. 但是，随着采样率的增加，ScGrassMC算法的收敛速度没有明显变化，这是由于ScGrassMC算法以子空间迭代法（subspace iteration）为基础，在Grassmann流形上使用非规范化度量，通过向目标函数的梯度引入缩放因子，使迭代过程能够快速收敛. 其余3种算法的收敛速度变化明显，其中IALM算法变化最为明显，在采样率为50%时，其收敛速度仅次于ScGrassMC算法.

实验2 　EDM恢复误差比较.

1）情况1下的EDM恢复误差. 如图2（a）所示为不同采样条件下EDM恢复误差的变化情况. 可以看出，在无噪声条件下，各算法的恢复误差均随采样率的增大而迅速减小，直到接近0. 当采样率为20%时，ScGrassMC算法的恢复误差接近0，算法性能优于另外3种算法；当采样率为30%时，LBD_MDS、ScGrassMC的恢复精度大致相同，优于IALM、OptSpace.

2）情况2下的EDM恢复误差. 各算法的野值噪声比（outlier ratio）均设为5%. 如图2（b）所示，LBD_MDS算法能较好地消除高斯噪声和野值噪声的影响，当采样率为30%时，能达到接近0的恢复误差；其余3种算法的补全精度受噪声影响较大，即使以90%的采样率进行采样，OptSpace算法的恢复误差也较大，不能有效处理高斯噪声和野值噪声.

3）情况3下的EDM恢复误差. 各算法的脉冲噪声比均设置为10%（以EDM行受到噪声污染为例），如图2（c）所示，与情况1相比，ScGrassMC的性能明显降低，而LBD_MDS和IALM的恢复误差仅有小幅增长. 当采样率超过30%时，LBD_MDS算法在脉冲噪声的影响下仍能够以较高的精度进行EDM补全. 因此，LBD_MDS具有较好的脉冲噪声容错性.

4）情况4下的EDM恢复误差. 野值噪声比设置为5%，脉冲噪声比设置为10%（以EDM行受到噪声污染为例）. 如图2（d）所示为算法在情况4下的EDM恢复误差，由于EDM受到复合噪声的污染，基于IALM、OptSpace和ScGrassMC的方法的恢复误差较大，不能对EDM进行较精确的估计. 与情况2、3相比，复合噪声条件下LBD_MDS算法的恢复误差有一定程度的增大，但在采样率为30%时仍具有较好的EDM补全性能.

实验3 　定位误差及其方差的比较. 分析信标节点的数量对节点定位误差的影响，将信标节点数量的变化范围设为3~20，采样率设置为50%. 如图3所示为在情况3、4下信标节点数量与定位误差之间的关系. 图中，Nb为信标节点的数量. 可以看出，随着信标节点数量的增加，各算法的定位误差均呈下降趋势，当信标节点数小于6时，定位误差随信标节点数量的变化较明显，因此，本研究将信标节点的数量设置为6. 另外可以看出，相比于另外3种算法，LBD_MDS算法的定位精度较高.

不同噪声条件下的定位误差及其方差随采样率的变化情况如图4、5所示. 与实验2对比可知，定位误差及其方差与EDM恢复误差的变化一致；相对于另外3种算法，当采样率大于15%时，LBD_MDS仍能以较低的误差和误差方差进行节点定位. 值得注意的是，当采样率大于30%时，LBD_MDS的定位误差及其方差不到另外3种算法的1/4. 当采样率为30%时，在情况2条件下，不同规模WSN中各算法性能的比较如表1所示. 表中，t为算法运行所需时间. 可以看出，各算法的定位误差随网络规模的增加而减小；在相同规模的WSN中，LBD_MDS具有较低的定位误差.

实验4 　定位误差累积分布比较. 如图6所示为当采样率为30%时，各算法在不同噪声条件下的节点定位误差累积分布. 如图6（a）所示，在情况2下，LBD_MDS定位误差小于1 m的概率高于95%，而IALM、OptSpace、ScGrassMC算法的相应概率均低于60%；如图6（b）所示，在情况3下，LBD_MDS定位误差小于1 m的概率为100%，IALM、ScGrassMC算法的相应概率均小于80%，而OptSpace的相应概率约为25%. 因此，LBD_MDS算法具有较好的定位性能.

实验5 　噪声污染程度对算法性能的影响比较. 为了评估EDM受噪声污染程度对各算法定位性能的影响，采样率固定为50%，野值噪声比的范围为0.05~0.50，脉冲噪声比的范围为0.1~0.5. 节点定位误差随噪声污染程度的变化趋势如图7所示. 图中，or为野值噪声比，pr为脉冲噪声比. 可以看出，随着噪声比的逐渐增加，各算法的定位误差均有不同程度的增大，相比于另外2种算法，LBD_MDS和IALM在不同的噪声污染程度下均具有较好的鲁棒性，其中LBD_MDS定位误差最小，对野值噪声和脉冲噪声的噪声容错能力较强.

提出基于Bregman散度矩阵补全的WSN定位算法，即LBD_MDS算法，用于解决脉冲噪声条件下WSN节点定位问题. 通过与IALM算法、OptSpace算法和ScGrassMC算法的对比实验可以看出，在不同噪声条件下，LBD_MDS算法在保证高效性的同时，在定位精度和鲁棒性方面优于其他算法，特别是当采样率达到一定程度（例如大于30%）时，LBD_MDS算法的定位误差不到另外3种算法定位误差的1/4，能够有效消除脉冲噪声对定位精度的影响. 然而，在无噪声条件下，LBD_MDS算法的定位精度不高；另外，与ScGrassMC算法相比，在算法收敛速度方面有待进一步提高.

在接下来的工作中，将进一步提高算法在无噪声条件下的定位精度，并在理论上对算法的收敛性进行严格分析，进一步提升算法的收敛速度. 另外，将尝试展开LBD_MDS算法在传感器网络数据收集、推荐系统和链路预测等领域的应用研究.