針灸療效評估中激光散斑襯比圖像穩(wěn)定性研究

仇福康，陳曉榮，劉宏業(yè)，衛(wèi)曉翔，張凱凱，孔 平

（1．上海理工大學(xué) 光電信息與計(jì)算機(jī)工程學(xué)院，上海 200093；2．上海健康醫(yī)學(xué)院 文理教學(xué)部, 上海 201318）

引 言

中醫(yī)理療如今被越來越多的人關(guān)注與接受。中醫(yī)理療的關(guān)鍵是通過外在物理刺激，如針灸、拔罐、推拿等，作用于人體氣血淤積之處，進(jìn)而使得人體的氣血得以通暢，疼痛得以緩解。中醫(yī)不同于以生物、物理、化學(xué)等為基礎(chǔ)的西醫(yī)，中醫(yī)的治療效果不能通過一些生理指標(biāo)來直接評估，一般是通過患者的反饋加之主治醫(yī)師的醫(yī)療經(jīng)驗(yàn)得出綜合評判。近期的研究表明，中醫(yī)理療的效果可以從皮膚表層的微循環(huán)的變化得以觀測[1]，而皮膚微循環(huán)的具體情況需要現(xiàn)代醫(yī)學(xué)儀器來探究。激光多普勒[2]是常用的一種應(yīng)用于監(jiān)測微循環(huán)系統(tǒng)的技術(shù)手段，現(xiàn)在市面上已經(jīng)存在開發(fā)成熟的激光多普勒血流儀。激光多普勒的作用原理是：在激光的單點(diǎn)照射下，生物的被測區(qū)域中流動的紅細(xì)胞會改變散射回去的激光的頻率，這些回到接收端的激光的頻率存在頻率偏移[3]，通過數(shù)學(xué)分析這些偏移量可以推導(dǎo)出被測血流的速度分布情況。激光多普勒具有侵入深度深、成像清晰的優(yōu)勢，但不足之處是單點(diǎn)掃描采集相應(yīng)的頻移信息，其成像速度過于緩慢[4]。

隨著研究的進(jìn)一步深入，檢測效果優(yōu)于激光多普勒單點(diǎn)掃描成像的激光散斑襯比成像[5](laser speckle contrast imaging, LSCI)的技術(shù)標(biāo)準(zhǔn)得以制定。相對于激光多普勒，LSCI具有成像速度快、成像維度為二維以及成像空間分辨率高等優(yōu)點(diǎn)[6]，這些成像優(yōu)勢使得LSCI可用來實(shí)時(shí)監(jiān)測針灸過程中被測對象皮膚表層的微循環(huán)。LSCI的原理是：在全場照明下，被測區(qū)域微循環(huán)中流動的紅細(xì)胞會產(chǎn)生散射激光的相移[7]，在一定曝光時(shí)間內(nèi)不同相移的激光在接收端的像平面產(chǎn)生不同程度的干涉，分析不同程度的干涉所形成的散射圖像即可得出被測微循環(huán)區(qū)域的血流速度分布。LSCI的分析是基于散斑圖像空間信息的數(shù)理統(tǒng)計(jì)，定義了反映血液中紅細(xì)胞流動速度信息的襯比值K[8]。在實(shí)驗(yàn)研究中發(fā)現(xiàn)，激光散斑襯比成像的襯比值K的穩(wěn)定性會隨著光源調(diào)整而變化，為了提高激光散斑襯比成像中流速信息反饋的準(zhǔn)確性，對影響散斑成像的襯比值及光源進(jìn)行測試與分析，并進(jìn)一步改進(jìn)散斑成像的襯比值的處理算法，以提高激光散斑成像在針灸效果評估中的穩(wěn)定性。

1 激光散斑襯比成像理論與處理算法改進(jìn)

1.1 成像理論基礎(chǔ)

LSCI的理論基礎(chǔ)是隨機(jī)波動的激光光強(qiáng)在CCD相機(jī)設(shè)定的曝光時(shí)間內(nèi)在探測器上形成的不同的光強(qiáng)積分，即成像面上對應(yīng)位置產(chǎn)生隨機(jī)分布的灰度值。在設(shè)定的曝光時(shí)間內(nèi)，激光光源保持穩(wěn)定輸出，成像面的灰度值分布反映了被測區(qū)域的運(yùn)動情況，靜止的部分顯示為黑白分明的散斑圖案，運(yùn)動部分則根據(jù)速度的快慢顯示為模糊程度不一的灰度圖案，該灰度圖案稱之為散斑圖[9]。定義散斑襯比K的數(shù)學(xué)表達(dá)式為

式中：β是一個(gè)常系數(shù)，其值由光強(qiáng)探測器探測單元的尺寸與散斑尺寸的比值決定；σI和表示空間滑動窗口內(nèi)的光強(qiáng)分布的標(biāo)準(zhǔn)差與均值；x是CCD相機(jī)曝光時(shí)間T與光強(qiáng)波動的自相關(guān)時(shí)間τ的比值，其中光強(qiáng)波動的自相關(guān)時(shí)間τ與被測區(qū)域內(nèi)散射介質(zhì)運(yùn)動速度v成反比。

1.2 處理算法改進(jìn)

激光散斑成像理論是基于散斑圖像灰度值的空間數(shù)理統(tǒng)計(jì)，散射激光的干涉分布的散斑襯比K可表示為

式中Ii為空間窗口內(nèi)某點(diǎn)的光強(qiáng)。散斑襯比K的定義是，在N×N的空間滑動窗口內(nèi)，激光光強(qiáng)的標(biāo)準(zhǔn)差與激光光強(qiáng)的均值的比值。圖1顯示了激光散斑空間算法的基本操作流程：滑動窗口（一般選取5×5大小）即在散斑圖像上自左至右掃描，分別計(jì)算掃描過的區(qū)域內(nèi)的圖像灰度值的標(biāo)準(zhǔn)差與均值的比值，分別記為K11、K12、K13、…；以此類推，掃描完散斑圖像的第一行，圖1中的滑動窗口下移一行，完成對散斑圖像第二行的掃描，所得襯比分別記為K21、K22、K23、…。

圖1 激光散斑空間算法的基本原理Fig.1 Basic principle of laser speckle space algorithm

由式(2)以及圖1所示的激光散斑空間算法，處理激光在CCD像面上形成的包含被測區(qū)域內(nèi)散射介質(zhì)的運(yùn)動信息的灰度值分布圖，即在一定曝光時(shí)間內(nèi)，運(yùn)動的散射介質(zhì)形成散射激光相移，進(jìn)而干涉形成光強(qiáng)的積分分布圖。基于空間散斑算法的基本公式，計(jì)算襯比分布的耗時(shí)隨著CCD采集的圖片像素的增大而增加。針灸療效的檢測系統(tǒng)需具備快速數(shù)據(jù)處理功能，Tom等發(fā)表的快速處理算法能有效縮短散斑運(yùn)算的耗時(shí)，記為Roll算法[10]，算法的基本表達(dá)式為

采用Roll算法能顯著減少運(yùn)算時(shí)間，例如采用 Intel? CoreTMi7-4790K CPU@4.00 GHz、內(nèi)存16 GB、Microsoft windows 7、操作系統(tǒng)為64位的計(jì)算機(jī)，在散斑圖像大小為1040×1392、滑動窗口設(shè)置為5×5的情況下，激光散斑基本算法的耗時(shí)為76.512 s，而Roll算法耗時(shí)為6.027 s，運(yùn)算效率明顯提升。因此在采用Roll算法的基礎(chǔ)上，可以實(shí)現(xiàn)針灸過程中皮膚表層血流微循環(huán)的實(shí)時(shí)監(jiān)測。研究發(fā)現(xiàn)，雖然運(yùn)算速度得到提升，但是由于激光光源的非理想性，以及激光散斑本身具有隨機(jī)波動的特性[11]，因此待測區(qū)域內(nèi)的靜止部分仍會隨著時(shí)間進(jìn)行隨機(jī)的、細(xì)微的波動。將激光散斑應(yīng)用到針灸治療效果的評估時(shí)，被測區(qū)域的相對靜止部分能得到一個(gè)比較穩(wěn)定的散斑襯比，使得針灸刺激區(qū)域的細(xì)微反應(yīng)可以被觀察、記錄與分析。

為使針灸療效檢測系統(tǒng)達(dá)到實(shí)時(shí)顯示、輸出穩(wěn)定的需要，在Roll算法的基礎(chǔ)上，添加了基于時(shí)間的數(shù)理統(tǒng)計(jì)分析，以期得到穩(wěn)定的襯比輸出。改進(jìn)的流程如圖2所示，記此算法為Roll with time-sampling，簡稱RT。采用與Roll算法相同的計(jì)算機(jī)配置進(jìn)行實(shí)驗(yàn)，在散斑圖像大小為1040×1392、滑動窗口設(shè)置為5×5的情況下，RT算法耗時(shí)為 6.192 s。相對于 Roll算法，RT算法的運(yùn)算效率存在細(xì)微的下降，這是因?yàn)檩d入多組數(shù)據(jù)進(jìn)行對比運(yùn)算會產(chǎn)生延時(shí)。

圖2 RT 提升散斑圖像襯比穩(wěn)定的流程Fig.2 RT enhances the stability of the speckle image contrast

RT的處理思路是，通過基于快速數(shù)據(jù)處理的Roll算法，加入5個(gè)抽樣時(shí)刻t1、t2、t3、t4、t5，對應(yīng)的襯比值為k1、k2、k3、k4、k5，抽樣間隔分別為 Δt1、Δt2、Δt3、Δt4，總間隔為ΔT=Δt1+Δt2+Δt3+Δt4。

針灸效果評估要求實(shí)時(shí)顯示，故設(shè)定Δt1=5 s、Δt2=10 s、Δt3=15 s、Δt4=25 s，這樣既可以得到穩(wěn)定化的襯比值，又可以滿足測量系統(tǒng)對實(shí)時(shí)性的要求。

2 針灸測量系統(tǒng)與系統(tǒng)參數(shù)優(yōu)化

2.1 針灸測量系統(tǒng)

基于LSCI成像系統(tǒng)的針灸測量系統(tǒng)主要由激光光源輸出、圖像采集、數(shù)據(jù)分析3個(gè)模塊構(gòu)成，如圖3所示，其中，區(qū)域1為選自背景部分的靜態(tài)對照區(qū)域，區(qū)域2為以毫針刺入的位置為中心的刺激區(qū)域。激光經(jīng)光束擴(kuò)展裝置照射到針灸作用的區(qū)域，圖中橢圓形區(qū)域?yàn)镃CD相機(jī)的有效視場。激光經(jīng)被測區(qū)域后，以不同角度散射并被成像器件捕獲，成像器件包括一定帶寬的濾波片以及面陣CCD相機(jī)。數(shù)據(jù)分析用LabVIEW與MATLAB混合編制的分析程序來完成。

圖3 基于激光散斑襯比成像的針灸測量系統(tǒng)Fig.3 Acupuncture measurement system based on laser speckle contrast imaging

2.2 針灸測量系統(tǒng)參數(shù)優(yōu)化

針灸測量系統(tǒng)所得襯比穩(wěn)定性的影響因素主要為測量系統(tǒng)參數(shù)以及散斑成像系統(tǒng)參數(shù)，其中：散斑成像系統(tǒng)參數(shù)包括散斑大小[12]、散斑空間處理的滑動窗口大小[13]；成像系統(tǒng)參數(shù)包括CCD的靶面大小[14]、激光光源的功率[15]。隨著實(shí)驗(yàn)的進(jìn)一步深入，激光的波長作為一個(gè)新的影響因素被列入考察范圍。

實(shí)驗(yàn)是基于國內(nèi)外認(rèn)可的空間散斑襯比分析方法，測量被測試者手部微循環(huán)在針灸刺激下的變化情況。測量系統(tǒng)激光發(fā)生器選用了4種，分別為He-Ne激光器(波長632.8 nm，功率30 mW)、半導(dǎo)體激光器1(波長785.0 nm，功率50 mW)、半導(dǎo)體激光器2(波長532.0 nm，功率300 mW)、半導(dǎo)體激光器3(波長671.0 nm，功率300 mW)，通過對比實(shí)驗(yàn)篩選出能產(chǎn)生最佳襯比的激光光源。激光發(fā)射的光束經(jīng)由兩片具有一定間隔的透鏡所組成的光束擴(kuò)散系統(tǒng)，發(fā)散形成激光光柱照射到被測區(qū)域，實(shí)驗(yàn)對象為平整A4紙。

實(shí)驗(yàn)分為4組進(jìn)行，通過變量控制及添加激光光強(qiáng)衰減器，將激光功率統(tǒng)一調(diào)制為30 mW。通過分析襯比的數(shù)學(xué)統(tǒng)計(jì)特性來評判襯比的穩(wěn)定性，其中，選取襯比的方差為穩(wěn)定性的判斷指標(biāo)，實(shí)驗(yàn)結(jié)果見表1。

表1 激光光源測試結(jié)果歸納Tab.1 Summary of laser source test results

通過分析對比各光源下襯比值的方差大小，可以得出：采用半導(dǎo)體激光器1可使襯比穩(wěn)定性最佳。因此在針灸效果評估的實(shí)驗(yàn)中，采用波長為785.0 nm的半導(dǎo)體激光器作為光源。

3 實(shí)驗(yàn)結(jié)果與討論

圖4是針灸實(shí)驗(yàn)激光散斑圖，實(shí)驗(yàn)環(huán)境溫度為26 ℃，測試者靜坐于實(shí)驗(yàn)臺前，右手放在被測區(qū)域內(nèi)，保持靜止?fàn)顟B(tài)。實(shí)驗(yàn)時(shí)長共22 min，分為針灸前2 min、針灸進(jìn)行中10 min、針灸后10 min。實(shí)驗(yàn)采取的抽樣間隔為1 min，每次采集5張圖片，圖片采集間隔依次分為5 s、10 s、15 s、25 s。

圖4 針灸實(shí)驗(yàn)激光散斑圖Fig.4 Original laser picture of acupuncture experiment

分別通過Roll算法及RT算法處理測得的激光散斑圖片，所得的襯比偽彩圖如圖5所示。圖5中，(a)表示Roll算法處理針灸前第1個(gè)抽樣時(shí)刻的原始散斑圖所得襯比偽彩圖，(d)表示RT算法處理針灸前第1個(gè)抽樣時(shí)刻的原始散斑圖所得襯比偽彩圖；(b)表示Roll算法處理針灸過程中第4個(gè)抽樣時(shí)刻的原始散斑圖所得襯比偽彩圖，(e)表示RT算法處理針灸過程中第4個(gè)抽樣時(shí)刻的原始散斑圖所得襯比偽彩圖；(c)表示Roll算法處理針灸后第5個(gè)抽樣時(shí)刻的原始散斑圖所得襯比偽彩圖，(f)表示RT算法處理針灸后第5個(gè)抽樣時(shí)刻的原始散斑圖所得襯比偽彩圖。偽彩圖區(qū)間設(shè)置為0(藍(lán)色)～70(紅色)。

圖5 Roll和 RT 算法處理所得襯比偽彩圖對比Fig.5 Comparison pseudo color pictures of Roll and RT algorithm

圖5中：標(biāo)記“2”的區(qū)域?yàn)獒樉拇倘肴梭w后人體皮膚表層血流微循環(huán)的觀察區(qū)域；標(biāo)記“1”的區(qū)域?yàn)殪o態(tài)參考點(diǎn)，用來評價(jià)算法處理原始散斑圖所得襯比圖分布的穩(wěn)定程度；經(jīng)過Roll算法處理所得的(a)、(b)、(c)3幅偽彩圖之間出現(xiàn)了比較明顯的整體的數(shù)據(jù)波動；相對而言，經(jīng)過RT算法處理所得的(d)、(e)、(f)3幅偽彩圖之間波動比較平穩(wěn)。區(qū)域2的變化情況能直觀表現(xiàn)針灸的治療效果，在RT算法的處理下，將圖像的整體抖動進(jìn)行了優(yōu)化，進(jìn)而可以比較準(zhǔn)確地反映出區(qū)域2的襯比值波動情況。

圖6給出了時(shí)長為22 min的手背針灸測試結(jié)果穩(wěn)定性對比圖，其中：“+”表示RT算法處理原始散斑圖像區(qū)域1數(shù)據(jù)后所得襯比均值隨時(shí)間的分布曲線；“□”表示Roll算法處理原始散斑圖像區(qū)域1數(shù)據(jù)后所得襯比均值隨時(shí)間的分布曲線；“◇”表示Roll算法處理原始散斑圖像區(qū)域2數(shù)據(jù)后所得襯比均值隨時(shí)間的分布曲線；“○”表示RT算法處理原始散斑圖像區(qū)域2數(shù)據(jù)后所得襯比均值隨時(shí)間的分布曲線。

理想狀態(tài)下，激光光源保持持續(xù)的穩(wěn)定輸出，測試系統(tǒng)所用的背景板是處于絕對靜止?fàn)顟B(tài)，此時(shí)，靜止背景部分的襯比值為1(實(shí)際上由于被測試人員的手部在測試過程中不能保持長時(shí)間的相對靜止，即存在肌肉的無意識舒張運(yùn)動，會使底端的背景板產(chǎn)生細(xì)微的抖動)。襯比的分布區(qū)間為0～1，數(shù)值越小表示對應(yīng)區(qū)域內(nèi)運(yùn)動的散射介質(zhì)運(yùn)動越快；同理，數(shù)值越大表示對應(yīng)區(qū)域的散射介質(zhì)運(yùn)動越慢。經(jīng)過RT算法優(yōu)化處理之后，由圖6可見：2 min之后(0～2 min為針灸前)針灸區(qū)域2開始呈現(xiàn)襯比減小的趨勢；3～12 min為針灸中，針灸區(qū)域2保持襯比減小趨勢并達(dá)到極小值；13～22 min為針灸后，針灸區(qū)域2出現(xiàn)襯比增加的趨勢，并最終回歸到靜止?fàn)顟B(tài)時(shí)的襯比分布值。相對而言，Roll算法處理區(qū)域2的原始數(shù)據(jù)后，襯比分布中存在明顯的波動，在針灸前第2個(gè)時(shí)間點(diǎn)，針灸中第4、6、10個(gè)時(shí)間點(diǎn)，針灸后的第4個(gè)時(shí)間點(diǎn)均產(chǎn)生異常襯比，這表明了Roll算法處理原始散斑數(shù)據(jù)所得襯比的穩(wěn)定性欠佳。參考靜態(tài)對照區(qū)域1的襯比分布，相對于Roll算法，RT算法處理所得襯比的波動程度明顯小很多，相對平穩(wěn)，襯比波動值為1.73%(襯比隨時(shí)間分布的標(biāo)準(zhǔn)差)，小于3.47%(激光散斑近期研究中定義的襯比波動穩(wěn)定臨界值[16])，因此相對于Roll算法，RT算法能在一定程度上改善襯比分布的穩(wěn)定性。

圖6 Roll和 RT 處理所得襯比隨時(shí)間分布Fig.6 Distribution of speckle contrast over time obtained by Roll and RT algorithm

4 結(jié) 論

通過光源對照組的選擇以及相應(yīng)的散斑襯比值K的方差分析比較可知，波長785.0 nm、功率50 mW的半導(dǎo)體激光器產(chǎn)生的散斑襯比穩(wěn)定性最優(yōu)。采用對散斑襯比穩(wěn)定性能最優(yōu)的近紅外激光光源進(jìn)行實(shí)驗(yàn)，對時(shí)長為22 min的實(shí)驗(yàn)數(shù)據(jù)分別進(jìn)行Roll算法和RT算法處理。通過對比分析Roll算法、RT算法處理靜態(tài)背景區(qū)域1和針灸區(qū)域2的原始散斑圖像數(shù)據(jù)所得的襯比隨時(shí)間分布狀況，可以得出：RT算法處理所得的襯比分布具有更佳的去散斑圖像抖動的效果。