楊 莉, 胡國兵*, 胡學(xué)龍, 吳珊珊

(1. 金陵科技學(xué)院電子信息工程學(xué)院, 南京 211169; 2. 揚州大學(xué)信息工程學(xué)院(人工智能學(xué)院), 江蘇 揚州 225127;3. 南京信息職業(yè)技術(shù)學(xué)院電子信息學(xué)院, 南京 210023)

在雷達、無線通信、認(rèn)知無線電及通信偵察等信號處理場合, 下不完全Gamma函數(shù)逆的求解有著廣泛的應(yīng)用,如通信參數(shù)的設(shè)定[1-2]、信號檢測[3]、參數(shù)識別[4-5]及系統(tǒng)模型設(shè)計[6]等．下不完全Gamma函數(shù)為非初等函數(shù), 一般很難得到精確的解析解．近年來,下不完全Gamma函數(shù)的逆的求解大多采用數(shù)值解算法[7], 但數(shù)值解算法的收斂性與速度受初始值選擇的影響較大且迭代運算的計算成本高,故硬件實現(xiàn)較困難．顯然, 在實際信號處理尤其是工程應(yīng)用中, 獲取下不完全Gamma函數(shù)逆的近似解析解,對于簡化計算和提高理論推導(dǎo)結(jié)果的可解釋性及其簡潔程度具有重要意義．目前, 有關(guān)下不完全Gamma函數(shù)逆解析求解的算法鮮為見諸報道．Dohler等[1]將下不完全Gamma函數(shù)表示成廣義合流超幾何函數(shù), 獲得了尺度參數(shù)接近于0時下不完全Gamma函數(shù)的逆的一種近似解析解, 但當(dāng)尺度參數(shù)大于0時算法無法工作． 本課題組在前期工作中[8]提出一種基于極值理論(extreme value theory, EVT)的下不完全Gamma函數(shù)逆的近似解析解, 并將其應(yīng)用于多天線空時分組碼(space-time block code, STBC)碼型識別閾值的確定．該方法計算精度較高, 但需要對2個歸一化常數(shù)進行計算, 故計算成本偏高．本文擬將下不完全Gamma函數(shù)逆的求取與獨立同分布Gamma類隨機變量最大值的Gumbel近似聯(lián)系, 將求逆運算轉(zhuǎn)化為對歸一化系數(shù)的求取,并利用尾部等價原理, 提出一種較為簡單的近似解析解．

1 問題的提出

在認(rèn)知無線電系統(tǒng)中, 為了區(qū)分多輸入多輸出-正交頻分復(fù)用(multiple input multiple outputorthogonal frequency division multiplexing, MIMO-OFDM)系統(tǒng)所采用的STBC碼型, Eldemerdash等[4]提出一種基于恒虛警準(zhǔn)則的識別算法．該算法以任意兩根接收天線信號的互相關(guān)函數(shù)為依據(jù), 將識別統(tǒng)計量定義為

(1)

式中Nr為接收天線的個數(shù),τ為延時量,Fij(τ)為接收天線之間的互相關(guān)函數(shù),Fc(τ)為各接收天線之間的互相關(guān)函數(shù)之和．由式(1)可知, 識別統(tǒng)計量即Fc(τ)的最大值近似服從自由度為2的中心卡方分布．利用給定的虛警概率pfa計算得到相應(yīng)的識別閾值λ, 并通過比較V與λ的大小實現(xiàn)對STBC的碼型識別．識別閾值λ可通過求解方程

1-pfa=[P(Nc,λ/2)]N+v

(2)

得到, 其中N為樣本長度,v為循環(huán)前綴的長度(一般設(shè)為N/4),Nc=Nr(Nr-1)/2,P(k,x)為歸一化下不完全Gamma函數(shù), 可表示為

(3)

2 本文算法

2.1 算法原理

由式(1)可知,文獻[4]中用于區(qū)分MIMO-OFDM系統(tǒng)空時分組碼碼型的檢驗統(tǒng)計量V漸進服從中心卡方分布,其分布函數(shù)為

F(x;Nc,2)=P(Nc,x/2)．

(4)

顯然,該分布函數(shù)是歸一化下不完全Gamma函數(shù), 其形狀參數(shù)為Nc．于是,式(2)可以寫為

(1-pfa)1/(N+v)=F(λ;Nc, 2)．

(5)

閾值λ可通過對下不完全Gamma函數(shù)求逆得到,即

λ=F-1[(1-pfa)1/(N+v)]．

(6)

據(jù)文獻[10]可知, 如果獨立同分布隨機變量序列服從中心卡方分布,則其最大值的極限分布為廣義極值分布的Ⅰ型分布, 即Gumbel分布．根據(jù)文獻[11], 對于n個底分布為F(x;Nc, 2)的獨立同分布隨機變量Xi,i=1,2,…,n, 假設(shè)其最大值Mn=max{X1,X2,…,Xn}, 則存在歸一化常數(shù)an∈R,bn>0,使得

(7)

式中Pr(x)為取概率運算,G(x)=exp(-exp(-x))為Gumbel函數(shù)．由公式

(8)

可計算得到式(7)中的歸一化常數(shù)．顯然,歸一化常數(shù)bn的求取就是通過對下不完全Gamma函數(shù)求逆得到．

根據(jù)式(6)和式(8),構(gòu)造如下關(guān)系式:

F-1[(1-pfa)1/(N+v)]≈F-1(1-1/m),

(9)

式中m為正實數(shù)．通常在算法設(shè)計中, 虛警概率pfa?1, 于是可將(1-pfa)1/(N+v)進行一階泰勒近似, 得

(10)

由于(N+v)/pfa?1,o(x2)表示式(10)的余項為x2的高階無窮小, 所以式(10)中的二階及以上分量可忽略不計．若取m=int[(N+v)/pfa], 則式(10)可轉(zhuǎn)化為

(1-pfa)1/(N+v)≈1-m-1．

(11)

于是, 式(6)中定義的閾值

λ≈F-1(1-1/m)≈bm,

(12)

式中bm表示當(dāng)n=m時的歸一化常數(shù)．顯然, 若F(x)具有初等函數(shù)的形式, 則其逆可取得解析解,相應(yīng)的歸一化系數(shù)an,bn也存在解析解．然而,此處F(x)是下不完全Gamma函數(shù), 其形式中包含無窮積分項,通常無法獲得解析解．于是,下不完全Gamma函數(shù)的求逆問題可轉(zhuǎn)化為尋找與下不完全Gamma函數(shù)尾部等價的分布函數(shù), 且該分布函數(shù)的歸一化常數(shù)具有解析解．研究表明, Gamma隨機變量的尾部與廣義Weibull類分布的尾部等價,且兩種分布均屬于Gumbel類函數(shù)的最大值吸引場[1], 歸一化常數(shù)b′m為

b′m=2{ln(m/Γ(Nc))+(Nc-1)lnBm+[(Nc-1)2lnBm-(Nc-1)2ln (Nc-1)+Nc-1]/Bm},

(13)

式中參數(shù)

Bm=ln(m/Γ(Nc))+(Nc-1)ln(Nc-1)．

(14)

進一步地, 根據(jù)參數(shù)Nc及m, 以解析解的形式給出式(6)的下不完全Gamma函數(shù)的逆

λevt=b′m,

(15)

式中λevt即為采用本文算法計算出的識別閾值．

綜上, 本文算法步驟如下:

1) 設(shè)定輸入?yún)?shù)．設(shè)定樣本長度N, 接收天線個數(shù)Nr, 循環(huán)前綴長度v, 虛警概率pfa等參數(shù);

2) 計算輔助參數(shù)．由Nc=Nr(Nr-1)/2及m=int[(N+v)/pfa]得到輔助參數(shù)值Nc和m, 代入式(14)可得Bm;

3) 輸出閾值計算結(jié)果．根據(jù)等式λevt≈F-1(1-1/m)≈b′m及式(13), 計算得到閾值λevt．

2.2 算法的誤差分析

本文算法誤差的主要來源為式(10)中泰勒一階近似產(chǎn)生的誤差和式(13)中歸一化常數(shù)bm的近似所引起的誤差．具體情形如下:

1) 對于式(10), 按照泰勒級數(shù)展開的誤差估計公式, 若令pfa=x,f(x)=(1-x)α,α=1/(N+v), 則其n階近似的皮亞諾余項存在如下不等式:

(16)

2) 對于式(13)定義的歸一化系數(shù),有

Fn(a′nx+b′n)-G(x)=h(x)(bn-b′n)/lnn+o(1/lnn),

(17)

式中h(x)為定義在正實數(shù)域上的有界函數(shù)[1]．由文獻[1]可知,式(17)的收斂速度階數(shù)為(ln(lnn)/lnn)2．當(dāng)n大于100時, 本文算法中所用歸一化常數(shù)的收斂階數(shù)約為0.11; 隨著n的增大,收斂階數(shù)逐漸趨于零,但收斂速度均趨于平緩．

綜上, 本文算法的適用范圍為虛警概率較小(小于0.01), 接收天線個數(shù)適中(一般小于7), 且樣本數(shù)較大的場合(一般要求樣本數(shù)大于100)．

3 性能仿真與分析

為了評估本文算法的性能, 將本文提出的近似解析解在不同參數(shù)條件下進行仿真分析, 并與文獻[8]算法及數(shù)值解法[7]進行性能及計算復(fù)雜度對比．仿真所用的計算平臺為DELL(型號: XPS13), 其CPU為英特爾雙核i7-8550U(1.8 GHz), 操作系統(tǒng)為64位Windows 10．仿真中,以相對誤差|λnum-λevt|/λnum為性能評估的指標(biāo), 其中λnum為數(shù)值解方法計算所得精確閾值,相對誤差越大則意味著相應(yīng)算法的準(zhǔn)確度越低．

3.1 算法的性能

圖1為當(dāng)樣本長度為1 024時,虛警概率分別取10-3,10-5,10-7, 在不同接收天線個數(shù)條件下, 分別利用本文提出的近似解析算法、數(shù)值解法及文獻[8]算法計算得到的STBC識別閾值如圖1所示．由圖1可見: 1) 當(dāng)虛警概率相對較小且接收天線個數(shù)適中時,由本文算法得到的近似解與數(shù)值精確解之間能較好地吻合, 相對誤差較小, 且本文算法的性能優(yōu)于文獻[8]算法; 2) 當(dāng)虛警概率變大時, 若接收天線個數(shù)較小 (不大于5), 由本文算法計算得到的識別閾值與數(shù)值解可較好吻合,相對誤差較小, 其性能優(yōu)于文獻[8]算法;若接收天線個數(shù)較大,本文算法的性能略次于文獻[8]算法; 3) 當(dāng)虛警概率恒定不變時,本文算法得到的近似解與數(shù)值精確解之間的相對誤差隨接收天線個數(shù)的增加而變大,其原因是式(13)給出的歸一化常數(shù)的近似解析解是通過尾部等價及有限項級數(shù)求和得到,當(dāng)下不完全Gamma函數(shù)的形狀參數(shù)(與接收天線個數(shù)呈正相關(guān))較大時,該近似式的誤差變大．

表1為假定信號樣本長度為1 024,接收天線個數(shù)為4, 虛警概率為0.001時,在相同條件下分別采用文獻[8]算法、數(shù)值解法[7]及本文算法對碼型識別閾值進行計算時的相對誤差及其運算時間統(tǒng)計結(jié)果．由表1可見:本文算法的相對誤差約為文獻[8]方法的1/3;運算時間約為文獻[8]方法的1/4,約為數(shù)值解法的1/240．這是由于本文算法僅須計算1個歸一化系數(shù),故耗時最少;而文獻[8]須計算2個歸一化系數(shù),運算量略有增加;數(shù)值解法則因需要進行多次迭代運算, 故復(fù)雜度最高．

表1 各算法的性能及運算時間

3.2 算法性能的影響因素

1) 虛警概率．當(dāng)樣本長度為1 024時,不同虛警概率條件下利用本文算法計算得到的識別閾值的相對誤差如表2所示．由表2可見, 當(dāng)接收天線個數(shù)保持不變時,本文提出的近似解析解的相對誤差隨虛警概率的增加而變大．其原因在于當(dāng)樣本長度、循環(huán)前綴長度以及接收天線個數(shù)等參數(shù)不變時,虛警概率越小, 等效m值越大, 本文算法中Gumbel函數(shù)的歸一化系數(shù)的收斂速度越快．

表2 虛警概率對碼型識別閾值相對誤差的影響

2) 信號樣本數(shù)．當(dāng)虛警概率pfa=10-4時, 本文算法的相對誤差與觀測樣本數(shù)之間的關(guān)系如圖2所示．由圖2可見: 1) 當(dāng)接收天線個數(shù)不變時,隨信號的觀測樣本數(shù)增加,本文算法的相對誤差變?。@是由于觀測樣本數(shù)變大導(dǎo)致等效m值變大, 本文算法中Gumbel函數(shù)的歸一化系數(shù)的收斂速度加快, 且尾部等價的近似性能越好,從而歸一化常數(shù)b′m的計算誤差相應(yīng)變小; 2) 當(dāng)樣本數(shù)固定時,本文算法的相對誤差隨天線個數(shù)的增大而變大．這是因為當(dāng)接收天線個數(shù)較大時歸一化常數(shù)b′m的近似誤差變大,從而導(dǎo)致算法的相對誤差變大．

4 結(jié)論

本文以MIMO-OFDM系統(tǒng)中空時分組碼碼型識別閾值計算為例, 研究了下不完全Gamma函數(shù)逆的求解問題,提出一種基于EVT及尾部等價理論的近似解析解．仿真結(jié)果表明: 當(dāng)虛警概率較小、天線個數(shù)適中且樣本量適中時,本文算法近似誤差較小,且計算復(fù)雜度優(yōu)于現(xiàn)有算法．本文算法便于在DSP/FPGA等硬件平臺實現(xiàn),對于提高算法工程實現(xiàn)的效率及結(jié)論的可解釋性具有重要意義,且該算法經(jīng)適當(dāng)修正后可用于其他涉及下不完全Gamma函數(shù)逆的應(yīng)用場合．后續(xù)將研究形狀參數(shù)較大時,下不完全Gamma函數(shù)逆的近似算法及其性能評估問題．