趙會夢，李臣鴻

（中南民族大學(xué)生物醫(yī)學(xué)工程學(xué)院膜離子通道與藥物研發(fā)實(shí)驗室，武漢 430070）

0 引言

隨著新型冠狀病毒肺炎疫情的暴發(fā)，醫(yī)護(hù)人員所必需的厚重防護(hù)服無疑會導(dǎo)致靜脈穿刺失敗率增加，這可能會延誤搶救患者的時機(jī)，甚至?xí)＜盎颊叩纳?，目前沒有特定的醫(yī)療輔助工具可以解決這一問題。因此，研究靜脈圖像的識別及其處理，定位血管的三維位置對于輔助靜脈穿刺等醫(yī)療手段具有重要意義。

血液中血紅蛋白對近紅外光的吸收特性使得近紅外光在血管相關(guān)的臨床應(yīng)用范疇中具有獨(dú)特優(yōu)勢，國內(nèi)外均積極探索近紅外光在血管顯影相關(guān)領(lǐng)域的應(yīng)用[1]。中國科學(xué)院西安光學(xué)精密機(jī)械研究所以其精湛的近紅外成像技術(shù)成功研制出國內(nèi)首臺高速投影式血管顯像儀，成像速率達(dá)25幀/s[2]。該研究對血管進(jìn)行了二維平面定位，但缺乏血管的深度信息[2]。拉脫維亞電子與計算機(jī)科學(xué)研究所的學(xué)者們構(gòu)建了用于血管結(jié)構(gòu)分析的近紅外成像系統(tǒng)，該系統(tǒng)能夠測試適合掌靜脈成像的近紅外波長，但是同樣缺乏血管的深度信息[3]。

提供二維信息的靜脈顯像儀的發(fā)展到目前為止已經(jīng)頗為成熟，能提供三維信息的設(shè)備卻仍然處在各大高校及研究機(jī)構(gòu)的研發(fā)之中，并沒有足以碾壓其他產(chǎn)品的優(yōu)勢產(chǎn)品。本研究希望為血管的三維定位研究提供新的方法和思路，進(jìn)而為醫(yī)護(hù)人員靜脈穿刺提供參考，達(dá)到輔助醫(yī)療的目的。因此，計算血管的皮下深度是本研究將靜脈信息從二維升級為三維的重點(diǎn)，通過探尋不同波長下靜脈圖像中靜脈灰度值和靜脈皮下深度之間的關(guān)系，建立模擬方程，以此實(shí)現(xiàn)提供血管的三維位置的目的。

1 近紅外靜脈圖像采集系統(tǒng)的搭建

近紅外靜脈圖像采集系統(tǒng)分為2個部分：第一部分為圖像采集設(shè)備。使用無畸變低照度近紅外定制攝像頭，加入700～1 100 nm的寬動態(tài)濾光片，能夠濾除可見光的干擾。第二部分為自行設(shè)計制作的多波長近紅外光源，該光源需要滿足的條件如下：（1）多波長近紅外光源在同一塊電路板上；（2）不同波長光源位置分布均勻，以免造成誤差；（3）不同波長的光源自動切換，以免手動切換導(dǎo)致誤差?；谝陨显O(shè)計考量，本研究中多波長近紅外光源由單片機(jī)控制，每隔3 s切換下一個波長，單片機(jī)型號為STM32F103c8T6。光源的燈珠購自深圳市超越光電有限公司，選取了780、805、825、850、880和940 nm 6種波長[4]以用于后續(xù)多波長近紅外圖像處理。多波長近紅外光源燈珠擺放位置示意圖及其實(shí)物圖如圖1所示，近紅外靜脈圖像采集系統(tǒng)實(shí)物圖如圖2所示。

圖1 多波長近紅外光源燈珠擺放位置示意圖及其實(shí)物圖

圖2 近紅外靜脈圖像采集系統(tǒng)實(shí)物圖

2 二維血管增強(qiáng)與提取

由于近紅外靜脈圖像不夠清晰，會影響到對血管邊界的提取，本研究先對850 nm波長近紅外光下拍攝的靜脈圖像進(jìn)行圖像增強(qiáng)[5]，以利于下一步將靜脈與背景進(jìn)行分割，明確靜脈的二維位置，然后細(xì)化出靜脈中心線，更好地輔助醫(yī)護(hù)人員進(jìn)行靜脈穿刺。本研究進(jìn)行圖像處理的平臺為MATLAB。

2.1 圖像增強(qiáng)

對比度受限自適應(yīng)直方圖均衡化（contrast limited adaptive histogram equalization，CLAHE）處理對于醫(yī)學(xué)領(lǐng)域的圖像，尤其是醫(yī)學(xué)紅外圖像的增強(qiáng)效果尤為明顯。本研究中CLAHE程序運(yùn)行時間（11～16 s）由tic、toc（MATLAB內(nèi)的計時程序，可以定量計算完成制定程序所消耗的時間資源）顯示。手背靜脈原圖如圖3（a）所示，經(jīng)過CLAHE處理后如3（b）所示。

2.2 血管與背景的分割

為了把血管輪廓完全從背景中分割出來，本研究對圖像進(jìn)行了二值化處理。二值化處理目的是將一定范圍內(nèi)的灰度突出成灰度級內(nèi)最大值，將其他所有灰度降低至0。二值化處理可以簡化圖像，減小圖像所包含的數(shù)據(jù)量，能夠很方便地突出感興趣目標(biāo)的輪廓[6]。本研究通過隨機(jī)采集感興趣區(qū)域的灰度值，取其平均值作為閾值來實(shí)現(xiàn)二值化。二值化處理后的近紅外靜脈圖像如圖3（c）所示。

觀察圖3（c）可以看出，血管的檢測情況良好，但是血管邊緣呈鋸齒狀，毛刺頗多，所以對該圖做了一些形態(tài)學(xué)處理[7]。形態(tài)學(xué)處理中，膨脹和腐蝕是常用的方法[8-9]。本研究使用5×5的矩形算子對二值化處理結(jié)果圖先進(jìn)行開運(yùn)算再進(jìn)行閉運(yùn)算，處理結(jié)果如圖3（d）所示。由結(jié)果圖比對可以看出，形態(tài)學(xué)處理后血管邊緣細(xì)小毛刺得到平滑，而且血管形態(tài)流暢，完整度高。

2.3 細(xì)化處理

圖像的細(xì)化通常是指將一幅二值化圖像骨架化的操作。本研究中對圖3（d）做細(xì)化處理，即提取血管的中心線條軌跡，將血管部分細(xì)化成只有一個像素寬度的線條，如圖3（e）所示。細(xì)化結(jié)果可以顯示血管走向，有助于標(biāo)識出靜脈穿刺部位。

圖3 二維處理前后的近紅外靜脈圖像

3 多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型的建立與計算

手背靜脈三維定位的關(guān)鍵在于血管深度的獲取，以下詳細(xì)介紹通過構(gòu)建多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型計算血管深度并進(jìn)行驗證的過程。

3.1 手背皮膚光學(xué)模型的建立

比爾-朗伯定律指出，電磁波穿過物體時被吸收的量與光程中存在的吸光分子數(shù)目成正比[10]。具體吸收量與被穿透物體的吸光系數(shù)、厚度等因素有關(guān)[11]。吸光系數(shù)即單位光程中存在的吸光分子數(shù)目的具象表示。比爾-朗伯定律適用于所有的電磁輻射和所有的吸光物質(zhì)，包括氣體、固體、液體、分子、原子和離子[12]。因此，本研究在開發(fā)血管深度計算算法時，參考了該定律，在此基礎(chǔ)上結(jié)合實(shí)際情況，建立了不同波長近紅外光下的手背皮膚光學(xué)模型。

皮膚位于人體表面，人身體各處皮膚形態(tài)頗有差異，本研究主要以手背靜脈網(wǎng)上方皮膚為例。本研究中手背皮膚光學(xué)模型的分層主要可以分成3層：上皮層、真皮層和皮下脂肪層，透過脂肪再向下就是本研究的待測目標(biāo)——手背靜脈。在該手背皮膚光學(xué)模型中，血管位于皮下的深度即上皮層、真皮層和皮下脂肪層厚度之和。由皮膚分層可以看到，該深度的波動主要由皮下脂肪層厚度的改變導(dǎo)致，所以本研究默認(rèn)上皮層和真皮層的厚度固定不變，皮下脂肪層厚度上下波動。

在同一波長的近紅外光垂直照射下，上皮層、真皮層和皮下脂肪層的吸收系數(shù)各有不同，在不同波長的近紅外光垂直照射下，各層吸收系數(shù)也不相同。近紅外光于皮膚各層的吸收系數(shù)是很多研究方向都需要用到的數(shù)據(jù)，所以在有規(guī)律可循的層面已經(jīng)有學(xué)者給出了計算方法，沒有規(guī)律的層面則通過測量得到。如福建省激光精密加工工程技術(shù)研究中心的研究人員在進(jìn)行老化皮膚光分布的蒙特卡洛模擬時，總結(jié)了上皮層、真皮層和血管吸收系數(shù)的計算方法[13]。2018年，李博等[14]在與肥胖相關(guān)的近紅外光譜研究中測量了脂肪層在近紅外波段的數(shù)據(jù)，計算中用到的氧合血紅蛋白、水等數(shù)據(jù)來源于上海化工數(shù)據(jù)共享服務(wù)平臺。

在各波段的近紅外靜脈圖像上取血管中間的灰度值，取其平均值得到數(shù)據(jù)如下：780 nm灰度值為130，805 nm灰度值為127，825 nm灰度值為121，850 nm灰度值為100，880 nm灰度值為113，900 nm灰度值為94。根據(jù)圖像的對比度和近紅外光在靜脈中的傳播特性，從中挑選805、850、940 nm 3個波長為計算模型的建立提供數(shù)據(jù)。

3.2 多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型的建立

本研究中多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型搭建步驟如下：（1）對同一個待測血管搭建不同波長的近紅外光照射環(huán)境就會產(chǎn)生一個參數(shù)矩陣，因為皮膚模型分為3層，所以從多波長近紅外圖像中挑選3個波長的圖像用來產(chǎn)生3×3的參數(shù)矩陣。（2）在每個波長的近紅外靜脈圖像中分別獲取血管中間點(diǎn)的灰度值，經(jīng)處理后與上述參數(shù)矩陣組成增廣矩陣。（3）對該增廣矩陣進(jìn)行求解，得到各層的深度計算值，將超聲血管成像設(shè)備提供的血管深度及血管直徑作為參考進(jìn)行輔助運(yùn)算和誤差修正。（4）仍然使用該待測血管，另取其他部位樣本點(diǎn)進(jìn)行驗證運(yùn)算，計算取樣部位的血管深度，同樣使用超聲圖像進(jìn)行驗證。多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型如下：

式中，n為灰度級最大值（本文為255）；m為灰度值；d1、d2、d3、d4分別表示上皮層、真皮層、皮下脂肪層的厚度和血管直徑，單位均為分別表示波長為λ時上皮層、真皮層、皮下脂肪層和血液的吸收系數(shù)。

多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型的超聲血管圖及灰度值取樣示意圖如圖4所示，其中圖4（a）是用來作為標(biāo)準(zhǔn)的超聲血管圖，由圖可得：d4=0.141 cm，因此建立的方程組如下：

在d1+d2+d3=0.15 cm的約束條件下解該增廣矩陣并進(jìn)行誤差修正，對結(jié)果保留2位小數(shù)，可得d1=0.01 cm、d2=0.09 cm、d3=0.05 cm。

3.3 深度恢復(fù)計算結(jié)果

在MATLAB中確定805 nm波長的灰度圖像左上角坐標(biāo)為（338，755）、右下角坐標(biāo)為（352，821），根據(jù)圖4（b）中紅色矩形標(biāo)記所示的矩形樣本讀取樣本灰度值，使用該值計算脂肪層厚度dcalculate。脂肪層厚度dcalculate與d1、d2的和組成血管深度模擬值dsimulation，使用血管深度模擬值dsimulation與其在原圖中的坐標(biāo)數(shù)據(jù)（即位置信息）繪制三維散點(diǎn)圖（如圖5所示）。計算血管深度模擬值dsimulation與超聲顯示血管參考深度的誤差，誤差數(shù)據(jù)柱狀圖如圖6所示，誤差數(shù)據(jù)三維散點(diǎn)圖如圖7所示。

圖4 多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型的灰度值取樣示意圖及參照標(biāo)準(zhǔn)（超聲血管圖）

圖5 血管位于皮下深度三維散點(diǎn)圖

圖6 誤差數(shù)據(jù)柱狀圖

從圖5中可以看出，整體結(jié)果分布較為集中，第338行至339行波動較大。從圖6中可以看出，有90.85%的數(shù)據(jù)誤差≤0.015 cm，有96.92%的數(shù)據(jù)誤差≤0.021 cm。從圖7中可以看出，在1 005個驗證點(diǎn)中，誤差較大的點(diǎn)多位于樣本邊界處，推測取樣時可能超出了血管邊界，導(dǎo)致該行附近數(shù)據(jù)誤差較大。去除第338行和339行樣本數(shù)據(jù)，使用（340，755）到（352，821）共計871個樣本點(diǎn)重新計算誤差，得到的誤差修正柱狀圖如圖8所示。

圖7 誤差數(shù)據(jù)三維散點(diǎn)圖

圖8 誤差修正柱狀圖

從圖8中可以看出，去除第338行和339行樣本數(shù)據(jù)后，誤差最大值已經(jīng)由0.041 cm降至0.033 cm。誤差修正后，有91.50%的數(shù)據(jù)誤差≤0.015 cm，有97.82%的數(shù)據(jù)誤差≤0.021 cm，較修正前有少許提升。由誤差分布規(guī)律以及計算原理可以得出，越靠近血管中心軌跡，計算誤差越小。

4 結(jié)語

本研究通過建立多波長近紅外靜脈深度恢復(fù)數(shù)學(xué)模型，恢復(fù)了二維場景中丟失的深度信息，從而對手背靜脈進(jìn)行了三維定位。經(jīng)871個樣本點(diǎn)的計算和驗證，有97.82%的把握將血管深度測量的誤差控制在0.021 cm以內(nèi)。本研究為現(xiàn)行市場上的二維投影式靜脈顯像儀升級為三維提供了參考[15]，將來可能為疫情下的醫(yī)護(hù)人員和患者提供更多幫助。

本研究仍存在不足之處，例如并未完全脫離醫(yī)護(hù)人員實(shí)現(xiàn)智能自動靜脈穿刺，尚需人為參與等。因此在后續(xù)研究中，我們會致力于開發(fā)自動計算并判斷最佳穿刺位置的醫(yī)療器械，為其配備相應(yīng)的硬件，使其能夠按照靜脈定位數(shù)據(jù)做到自動配藥[16]、自動靜脈穿刺以及自動給藥，為醫(yī)護(hù)人員提供更大的幫助。