羅春海，甘官壽，黃子寧，鄭宜榮，黃開德，楊志勇

（1.佛山科學(xué)技術(shù)學(xué)院 數(shù)學(xué)與大數(shù)據(jù)學(xué)院，廣東 佛山 528000；2.南昌航空大學(xué) 軟件學(xué)院，江西 南昌 330000）

0 引 言

設(shè)備免持定位技術(shù)（Device-Free Localization, DFL）[1]最早由Moussa和Youssef在2007年提出，這是一種基于無線傳感網(wǎng)絡(luò)中射頻節(jié)點接收信號強度（Received Signal Strength, RSS）值的變化進行目標(biāo)定位的方法。2009年由Wilson和Patwari首次將DFL技術(shù)問題轉(zhuǎn)化為射頻層析成像（Radio Tomographic Imaging, RTI）[2]技術(shù)問題，開辟了新的的研究方向。

RTI技術(shù)需要在目標(biāo)環(huán)境周圍搭建無線傳感網(wǎng)絡(luò)（Wireless Sensor Networks, WSN）[3]，形成覆蓋目標(biāo)環(huán)境的監(jiān)測區(qū)域，由于該區(qū)域內(nèi)的物體會吸收、反射、繞射或散射部分傳輸能量，引起節(jié)點接收端的接收信號強度（Received Signal Strength, RSS）值變化，利用該變化對其進行陰影衰落分離，便可實現(xiàn)對興趣目標(biāo)成像的位置估計。因射頻信號具有可穿透非金屬物質(zhì)（例如墻體、煙霧及樹木等障礙物）、不受非微波（例如光線、聲音、溫度等）等因素干擾的優(yōu)越性，所以RTI技術(shù)在智能監(jiān)控、環(huán)境監(jiān)測、資源勘探、搶險救援以及軍事反恐等領(lǐng)域有著巨大的應(yīng)用前景。

在實際環(huán)境中，射頻鏈路RSS基準(zhǔn)值易受環(huán)境因素影響，興趣目標(biāo)與環(huán)境發(fā)生耦合效應(yīng)引發(fā)多徑干擾噪聲，導(dǎo)致陰影衰落分離困難，成像受到干擾。因此，建立合理的噪聲分布模型是目標(biāo)重構(gòu)的重要基礎(chǔ)。

本文提出一種基于魯棒學(xué)習(xí)與時空信息素抗干擾射頻層析成像方法。通過矩陣的低秩性完成對測量矩陣的混合高斯噪聲分解，進行迭代以得到背景鏈路矩陣，利用EM算法求解模型參數(shù)獲得前景信息。在此基礎(chǔ)上，根據(jù)目標(biāo)運動具有時空相關(guān)性的特點，引入時空信息素來有效抑制多徑干擾噪聲，提高目標(biāo)定位的精度和魯棒性。

1 基于魯棒子空間與時空信息素跟蹤方法

1.1 時序RSS數(shù)據(jù)建模

射頻層析成像應(yīng)用需要在目標(biāo)區(qū)域周圍部署射頻信號收發(fā)陣列，形成無線傳感網(wǎng)絡(luò)。其中，陣列由K個射頻信號收發(fā)節(jié)點組成，任意一對節(jié)點間的通信被視為一條射頻鏈路，共M=K(K-1)/2條射頻鏈路。

在WSN中，興趣目標(biāo)與其射頻環(huán)境會發(fā)生耦合效應(yīng)，該區(qū)域內(nèi)的物體（包含障礙物）會吸收、反射、繞射或散射鏈路的傳輸能量。部署射頻節(jié)點信號收發(fā)陣列如圖1所示。

圖1 部署射頻節(jié)點信號收發(fā)陣列

針對i條鏈路的接收信號強度RSS值，建立鏈路損耗模型如下：

其中：Pi為射頻節(jié)點發(fā)射功率，單位為dB；Li為射頻節(jié)點間距離、天線工作模式、設(shè)備阻抗等方面的靜態(tài)能量損耗，單位為dB；Si為射頻信號受物體阻擋產(chǎn)生的陰影衰減，單位為dB；Fi為物體在網(wǎng)絡(luò)中由多徑效應(yīng)引起的陰影衰落；εi為測量噪聲。

根據(jù)式（1），對t時刻的第i條鏈路RSS值yi建立數(shù)學(xué)模型：

其中：Pi(t)、Li(t)、SiLoss(t)、FiLoss(t)分別為t時刻的傳輸功率損耗、靜態(tài)能量損耗、障礙物與環(huán)境等因素引起的陰影衰減、障礙物與環(huán)境等因素引起的多徑衰落；SiTarget(t)為t時刻的興趣目標(biāo)相關(guān)陰影衰減；FiNoise、εi(t)分別為t時刻的動態(tài)環(huán)境因素引起的多徑干擾噪聲與測量噪聲。其中FiNoise隨時間變化且在噪聲中起主導(dǎo)作用，因此本文針對該種噪聲集合建立混合高斯噪聲模型。

在以往基于RSS衰減特征的RTI定位算法中，重構(gòu)圖像需要獲取RSS值變化鏈，針對射頻鏈路RSS基準(zhǔn)值易受環(huán)境因素影響的情況，本文構(gòu)造一個由RSS值yi按時間間隔Δt組成時序測量矩陣Y，表達式為：

其中：Δt=tj+1-tj，tj={tj|j=0, 1, ...,L}表示信號接收時刻；yi(tj)代表第i條鏈路在tj時刻的RSS值。因此對測量矩陣Y建模如下：

其中：Y∈RM×L為時序測量矩陣；U∈RM×r和VT∈Rr×L代表時序背景鏈路矩陣；X∈RM×L代表時序前景信息矩陣[4]。

tj時刻包含噪聲的前景信息矩陣為：

其中，xi(tj)為前景信息矩陣X在tj時刻的第i行元素，代表tj時刻第i條鏈路的前景信息，且xi(tj)為像素向量，等價于pi(tj)=(p1(tj)...pk(tj)...pd(tj))1×d，k為向量索引，d代表像素個數(shù)。將行向量pi(tj)擴展為像素矩陣P(tj)，可代表tj時刻目標(biāo)重構(gòu)影像分布情況，表達式為：

其中，pik(tj)代表tj時刻重構(gòu)圖像的第i條鏈路上任意一點k的像素取值。

利用低秩矩陣分解[5,6]將時序前景信息矩陣X映射到低維空間：

式中：⊙為Hadamard乘積；U表示基矩陣，U∈RM×r；V表示系數(shù)矩陣，V∈RL×r；r<

式中：ui、vj分別為U、V矩陣的第i個、第j個行向量，代表背景鏈路；xij代表前景信息，其中包括前景像素信息與噪聲信息；Φ代表相關(guān)陰影衰落模型矩陣，本文采用橢圓權(quán)重模型。

為盡可能地重構(gòu)前景信息，利用混合高斯（Mixture of Gaussian, MoG）[7]模型可擬合任何連續(xù)分布噪聲的特性來對前景信息進行建模，則前景信息矩陣X的稀疏信號可建模如下：

式中：xij為前景信息含噪聲合集；K是混合高斯個數(shù)，k表示其中任意一個分量編號；πk≥0為混合高斯權(quán)重，并且表示均值為0、方差為σ2的高斯分布；為第k個高斯分量的分布。

因此，測量矩陣Y的每個元素yij的概率可建模為：

式中：Π={πi|i=1, 2, ...,πK}；Σ={σi|i=1, 2, ...,σK}。

因此，測量矩陣Y的最大似然函數(shù)可寫為：

式中：Π為所有高斯分量權(quán)重集合；Σ為所有高斯分量方差集合；Ω為測量矩陣Y未丟失數(shù)據(jù)集合；i和j為集合Ω內(nèi)的行和列索引，之后的研究都在集合Ω內(nèi)展開。

1.2 基于魯棒子空間學(xué)習(xí)模型的參數(shù)估計

利用EM[8-9]算法求解有關(guān)U、V、Π、Σ參數(shù)的最大可能性解，則有：

（1）E步驟

定義測量數(shù)據(jù)相關(guān)隱變量：z→Y，令隱分布q(zijk)=γi,j,k表示混合高斯分布中每個數(shù)據(jù)來源于第k∈{1, 2, ...,K}個分布的概率，則隱變量有離散取值集合Z={z1,z2, ...,zK}，采用E步驟對隱變量yijk進行計算，即：

（2）M步驟

通過E步驟得到隱變量，后驗計算模型參數(shù)，M步驟計算框架可優(yōu)化為：

對Π、Σ參數(shù)進行迭代計算：

對U、V參數(shù)進行迭代計算：

根據(jù)低秩矩陣分解形式，權(quán)重矩陣W的元素wij為：

式（16）的極大化完全等價于加權(quán)L2LRMF問題，因此本文選擇實現(xiàn)簡單且性能良好的ALS（Alternated Least Squares）算法對U、V進行交替固定值求解，收斂的U、V矩陣為最優(yōu)解結(jié)果：

對測量矩陣Y使用MoG-LRMF模型，基于魯棒子空間學(xué)習(xí)模型的參數(shù)估計算法流程如下：

（1）輸入：測量矩陣Y，Y的非缺失項的索引集Ω。

（2）隨機初始化Π、Σ、U、V，MoG-LRMF模型參數(shù)K以及迭代閾值。

（3）重復(fù)。

（4）采用E步驟計算γijk(i=1, 2, ...,d；j=1, 2, ...,n；k=1,2, ...,K）。

（5）針對Π、Σ采用M步驟計算πk、σk2(k=1, 2, ...,K）。

（6）針對U、V采用M步驟利用ALS算法計算，其中W由式（17）計算得出。

（7）自動K調(diào)諧：如果存在i、j使得，則通過，將第i和第j個高斯分量合并為唯一的高斯分量，式中ni表示第i個高斯分量中的元素數(shù)，然后從Π、Σ移除第j個高斯參數(shù)。最后令K=K-1。

（8）直到收斂。

（9）輸出：U、V背景鏈路矩陣。

2 時空信息素抗干擾成像方法

根據(jù)以上方法得到的時序前景信息矩陣X重構(gòu)圖像，需要用到相關(guān)陰影衰落模型來刻畫像素值的分布情況，本文采用橢圓權(quán)重矩陣Φ來決定像素值的非零權(quán)值：

其中：Φik(tj)表示第tj時刻重構(gòu)像素權(quán)值，即tj時刻在第i條鏈路第k個像素pik(tj)的加權(quán)值；d為收發(fā)射頻節(jié)點間距；dik

Send(tj)和dikReceive(tj)為像素點pik(tj)與兩個收發(fā)節(jié)點的距離；λ為橢圓超出鏈路距離的長度，為調(diào)節(jié)橢圓寬度大小參數(shù)，通常設(shè)置為一個較小的值，且對所有測量鏈路所考慮的尺度大小保持一致。第i條鏈路橢圓權(quán)重模型如圖2所示。

調(diào)查結(jié)果表明，濕式誘捕器投放間距10、15、20 m處理誘捕茶尺蠖成蟲總量分別為52、316和100頭，平均每臺為13、79和25頭，可見間距15 m處理誘捕量最多，并與其他2個處理均達到顯著差異，較好地控制了茶園茶尺蠖。因此，濕式誘捕器在茶園的設(shè)置密度以間距15 m較好。

圖2 第i條鏈路橢圓權(quán)重模型示意圖

根據(jù)相關(guān)陰影衰落橢圓權(quán)重模型[10-11]，可以大致刻畫出某一時刻任意一個像素點與所在鏈路之間的關(guān)系，其中加權(quán)值則代表像素點處陰影衰落影響鏈路的程度，便可獲得重構(gòu)圖像：

其中：s(tj)=(tj)?h(tj)為tj時刻帶有時空信息素的重構(gòu)像素；h(t)為高斯卷積核，決定時空信息素點值的擴散程度；ρ為揮發(fā)率，代表時空信息素揮發(fā)程度，取值區(qū)間為[0，1]，其中ρ越大揮發(fā)速率越快，ρ越小揮發(fā)速率越慢。

3 實驗與結(jié)果分析

3.1 射頻感知網(wǎng)絡(luò)實驗環(huán)境搭建

本文的WSN實驗環(huán)境在室內(nèi)與室外均有部署，覆蓋監(jiān)測區(qū)域大小為5 m×5 m，一次性共部署20個射頻信號收發(fā)節(jié)點，相鄰節(jié)點間距為1 m。室內(nèi)室外實驗場景與平面圖如圖3所示。

圖3 室內(nèi)室外實驗場景與平面圖

網(wǎng)絡(luò)基于MICAz節(jié)點[12-13]，工作頻率為2.4 GHz，采用IEEE 802.15.4標(biāo)準(zhǔn)通信協(xié)議。節(jié)點的傳輸功率為-10 dBm。選擇1號節(jié)點作為基站，對其他射頻節(jié)點采集的數(shù)據(jù)通過USB接口傳輸至計算機后進行處理。

根據(jù)部署順序，對射頻節(jié)點逐個標(biāo)記ID號并設(shè)置相應(yīng)傳輸協(xié)議。為避免各節(jié)點間傳輸沖突，具體傳輸協(xié)議內(nèi)容如下：從節(jié)點ID為1至20發(fā)送數(shù)據(jù)，當(dāng)任意一個節(jié)點為發(fā)送信號狀態(tài)時，其余射頻節(jié)點切換為接收狀態(tài)，并檢查該發(fā)送節(jié)點ID編號；隨后按ID編號順序輪流切換收發(fā)狀態(tài)；重復(fù)檢查是否輪到當(dāng)前節(jié)點為發(fā)送狀態(tài)，否則節(jié)點為接收狀態(tài)并繼續(xù)等待。若出現(xiàn)信號阻塞、射頻節(jié)點響應(yīng)超時等情況，則立即跳過當(dāng)前節(jié)點，按照ID編號順序輪至相鄰節(jié)點繼續(xù)發(fā)送信號。

基站接收節(jié)點所收到的數(shù)據(jù)，包括發(fā)送節(jié)點的編號以及編號對應(yīng)節(jié)點的RSS值，若接收到基站本身的信號則RSS值設(shè)定為-45 dBm。一個節(jié)點傳輸一次數(shù)據(jù)的過程記為1個周期（1個周期花費時間約為0.545 s）。

3.2 實現(xiàn)過程

所有方法都在MATLAB 2020a上實現(xiàn)，并配備在Intel（R） i7-9750H（CPU）、GeForce GTX 1650（GPU）和16 GB RAM的PC上運行。

如圖4所示，本文射頻層析成像過程主要分為三個階段：（1）通過部署WSN節(jié)點，覆蓋目標(biāo)空間形成監(jiān)測區(qū)域，采集來自各個射頻鏈路的節(jié)點RSS值組成測量矩陣；（2）對測量矩陣進行低秩矩陣分解，利用混合高斯模型對包含多徑衰落噪聲的前景部分進行建模，利用魯棒子空間學(xué)習(xí)方法求解參數(shù)，實現(xiàn)對測量矩陣Y的有效分解；（3）對分解出的前景信息矩陣X，利用目標(biāo)運動具有時空相關(guān)性的特點，使用橢圓權(quán)重模型對鏈路上像素點進行加權(quán)表示，引入時空信息素圖抑制多徑干擾噪聲，提高目標(biāo)定位的精度和魯棒性，最后得到具有噪聲抗干擾性的重構(gòu)圖像。

圖4 基于魯棒學(xué)習(xí)與時空信息素抗干擾方法框架

3.3 結(jié)果分析

表1的重構(gòu)圖像結(jié)果表明，本文方法在室內(nèi)（外）針對單目標(biāo)運動的射頻層析成像結(jié)果清晰準(zhǔn)確，能夠呈現(xiàn)出基本貼近實際的目標(biāo)運動軌跡，并有效對抗多徑效應(yīng)干擾；在室內(nèi)復(fù)雜場景下，針對雙目標(biāo)同向走圓、反向走圓的復(fù)雜運動，依然能較為清晰地重構(gòu)出大致結(jié)果；本文方法在室內(nèi)（外）針對雙目標(biāo)運動的射頻層析成像結(jié)果較為清晰，并且在運動物體較為貼近時，仍能有效對抗物體運動間產(chǎn)生的多徑效應(yīng)干擾，在實際實時成像中呈現(xiàn)出分辨度更高的結(jié)果。

表1 本方法室內(nèi)（外）的單（雙）運動目標(biāo)重構(gòu)圖像結(jié)果

3.3.2 平均定位誤差分析

如圖5的重構(gòu)軌跡對比結(jié)果表明，本文方法在室內(nèi)（外）情況下，針對單（雙）運動目標(biāo)可以很好地重構(gòu)出貼合實際運動的目標(biāo)軌跡，重構(gòu)與真實軌跡誤差均在1 m范圍以內(nèi)。本方法能夠有效對抗多徑效應(yīng)干擾，大幅提高誤差精度，縮小定位誤差范圍，具有很好的抗干擾性與魯棒性。

圖5 室內(nèi)（外）的單（雙）運動目標(biāo)的實際運動軌跡與本方法的重構(gòu)軌跡對比結(jié)果

表2的實驗結(jié)果表明，在室內(nèi)（外）WSN覆蓋區(qū)域大小為5 m×5 m的監(jiān)測范圍內(nèi)，單運動目標(biāo)定位誤差最大不超過0.427 m，雙目標(biāo)運動定位誤差最大不超過0.640 m，整體平均定位誤差最大不超過0.162 m，最優(yōu)情況下單（雙）目標(biāo)平均定位誤差最小值僅為0.002 m。

表2 室內(nèi)（外）的單（雙）運動目標(biāo)平均定位誤差結(jié)果

圖6所示的室內(nèi)（外）單（雙）運動目標(biāo)平均定位誤差累計分布計算中采用真實位置與估計位置的OMAT距離來衡量定位性能，其計算公式為：

圖6 室內(nèi)（外）的單（雙）運動目標(biāo)平均定位誤差累計分布曲線

其中：g為真實目標(biāo)的位置；G為真實目標(biāo)的個數(shù)；Ψ表示重構(gòu)出的所有非零位置；OMAT 距離為從所有非零位置中選出最接近真實目標(biāo)位置的G個位置作為估計位置后計算與真實目標(biāo)位置g的誤差距離。

實驗結(jié)果表明，室內(nèi)（外）的單（雙）運動目標(biāo)在0.2 m處誤差均快速收斂，接近90%的概率使定位誤差不超過0.4 m。

3.3.3平均定位時間分析

表3的實驗結(jié)果表明，本方法針對室內(nèi)（外）的單（雙）運動目標(biāo)平均定位時間均在0.03～0.05 s左右，最快平均定位時間在0.0317 s，滿足毫秒級定位，具有快速性。

表3 室內(nèi)（外）的單（雙）運動目標(biāo)平均定位時間結(jié)果

4 結(jié) 語

針對射頻鏈路接收信號強度（RSS）基準(zhǔn)值易受環(huán)境因素影響，以及興趣目標(biāo)與環(huán)境耦合會引發(fā)多徑干擾造成射頻層析成像性能衰退問題，本文提出基于魯棒學(xué)習(xí)和時空信息素的抗干擾成像方法。鏈路RSS表示為環(huán)境相關(guān)的背景和目標(biāo)運動相關(guān)的前景兩部分，構(gòu)造RSS時間序列矩陣，其中背景部分通過矩陣低秩性建模，前景通過混合高斯分布建模，形成對RSS測量噪聲魯棒的MoG-LRMF模型，設(shè)計EM算法求解模型參數(shù)以獲得前景信息。以此為基礎(chǔ)實現(xiàn)目標(biāo)成像，根據(jù)目標(biāo)運動具有時空相關(guān)性的特點，引入時空信息素圖抑制多徑干擾噪聲，提高目標(biāo)定位的精度和魯棒性。實驗結(jié)果表明，本方法在室內(nèi)及室外均表現(xiàn)良好，對于存在外部物體運動或阻擋的情況下仍能有效實現(xiàn)單（雙）目標(biāo)的成像，具有較好的抗干擾性、魯棒性和快速性。