基于蒙特卡羅的腦血氧檢測光源—探測器最優距離研究

曹 焱，邢志明，趙 斌，董祥美，史以玨，高秀敏

（1.上海理工大學 光電信息與計算機工程學院，上海 200093；2.上海交通大學醫學院附屬瑞金醫院，上海 200025）

0 引言

大腦是人體的一個重要器官，像一個中樞系統一樣指揮著人們的日常行為與活動。腦部的血流量約占心臟血液輸出量的15%，耗氧量約占全身耗氧量的20%［1］，所以大腦對于缺血缺氧是十分敏感的［2］。臨床上對腦血氧的監測可以幫助醫生在手術中實時監測患者健康狀況［3-4］，以及進行一些腦部疾病的診斷［5-6］。因此，用于臨床監測和診斷的血氧監測儀器可為醫生提供很大的幫助。

利用功能性近紅外光譜技術［7］可以很容易地檢測到腦功能的血流動力學信息［8］。在實際的腦血氧測量中，被檢測者往往會因人體的生理活動（如心動周期、呼吸、血壓變化等）引入干擾，從而給腦血氧信號檢測造成了一定影響，這些生理干擾可以采用自適應濾波法［9-10］去除。自適應濾波法對于抑制特定的生理干擾有比較好的效果，但是需要借助額外的設備獲取參考信號。通常使用單個探測器進行腦血氧信號檢測時會面臨以上問題，因為僅有一個探測器無法獲取參考信號。為解決單個探測器檢測所面臨的問題，有學者提出了多距離測量方法，即使用雙探測器測量腦血氧信號［11-12］。一個探測器檢測深層組織中帶有生理活動干擾的血氧信號，另一個探測器檢測淺層信號作為參考信號，對兩個探測器檢測的信號作差分處理［13］，獲得純度更高的局部信號用于后續計算處理。此外，臨床和科研中光源—探測器（Source-Detector，SD）距離的選擇大多根據經驗公式［14］，即探測器能檢測到的最大組織深度是光源與探測器距離的一半［15］，或者參考市面上現有的腦血氧檢測產品來選擇光源—探測器之間的距離［16］。這些方法選擇的光源與探測器距離是固定的，但由于顱腦各部分組織結構針對不同的年齡、性別、個體等會存在差異，如正常人和頭部有損傷的病人（如硬膜腦血腫［17］）等，經驗公式可能無法很好地適用于所有人群的腦血氧檢測，無法滿足光源—探測器之間的最佳距離選擇需求。

為解決上述腦血氧檢測中SD 距離選取時遇到的問題，本文根據多距離測量方法要求，給光源—探測器選擇一組恰當的距離。光源放置于固定位置時，探測器需一個放置于近處，用于檢測頭皮和顱骨處的血液信號作為參考信號；另一個放置于遠處，用于檢測灰質處受到生理活動干擾的血液信號。再結合蒙特卡羅方法，從光在組織中的傳輸特性著手，同時考慮到顱腦組織的結構差異，選定一組顱腦組織的特性參數，建立5 層顱腦組織平板模型進行仿真。通過統計大量光子在組織中的傳輸情況，并對仿真得到的數據進行分析與處理，可以很好地克服在顱腦組織存在差異時應用經驗公式的缺陷，從而選擇出符合該組顱腦組織參數的最優的光源—探測器距離。

1 蒙特卡羅模型與方法

1.1 蒙特卡羅原理

生物組織在光學上可近似看作混沌介質，當光照射在生物組織上時，會發生反射、吸收和散射等現象，光在生物組織中的宏觀傳播行為已不能用基本的朗伯—比爾定律來描述。近些年得益于計算機技術的快速發展，蒙特卡羅（Monte Carlo）方法在各個物理仿真實驗中使用得越來越多，尤其在生物醫學、光學等領域［18］。Monte Carlo 模擬［19-20］是一種基于概率統計的數學物理模擬方法，通過相應的概率模型或隨機過程計算，模擬粒子在組織中傳輸時受到的影響［21］。通過對大量光子在組織中傳播路徑的仿真，追蹤光子運動的軌跡點，可得到光在組織內部和表面的統計分布規律。

1.2 光在組織中傳輸的蒙特卡羅模擬方法

本文利用多層蒙特卡羅（Monte Carlo Multi-Layer，MCML）程序進行仿真，MCML 程序是Wang 等［19］用標準C 語言開發的程序包，用于模擬分層組織的光子傳輸［22］。蒙特卡羅程序的基本流程如下［23-24］：

（1）確定起始跟蹤點和初始化光子（權重、方向、位置）后，在組織表面發射光子，并為其分配隨機步長。

（2）移動步長s 后，更新光子位置，確定光子下一次的碰撞位置。

（3）光子移動步長后會隨機進行吸收和散射，散射作用會改變光子的前進方向，吸收作用會降低光子權值。根據當前位置計算由于散射造成的方向偏轉，更新傳播方向余弦。計算光子在當前位置由于吸收造成的權重衰減，更新權重值，并對光子權重是否低于閾值進行判斷。如果低于閾值，或逸出組織上下表面時，采用輪盤賭來決定光子存活。判定為死亡的，結束當前光子的追蹤。

（4）對光子是否碰撞邊界進行判斷，根據判斷結果進行相應的處理。

（5）對于未逃逸出組織的光子繼續進行步驟（2）—（4），直到程序結束。

蒙特卡羅模擬流程如圖1所示。

Fig.1 Flow of Monte Carlo simulation圖1 蒙特卡羅模擬流程

1.3 光在腦組織中傳輸路徑的仿真研究

研究光在顱腦組織中的傳輸路徑，需對每個光子的傳輸軌跡進行追蹤，但Wang 原有的程序并不能實現此功能。因此，需要在Wang 的仿真程序中進行優化，以獲得光子在組織中的傳輸路徑。由于光子在組織內傳輸是隨機、無序的，所以采用數據鏈表來保存光子的傳輸路徑。光子在組織中傳輸時的平均自由路程很短，在運算時保存的信息會占用很大的計算機運行內存，故用d來判斷光子運動距離，優化保存的信息。

式中，d表示光子傳輸時相鄰兩個位置的距離，（xk，yk，zk）表示光子當前的位置，更新后的位置坐標為（x，y，z）。只有當d值大于某一閾值時，才將當前位置點的坐標記錄到鏈表中［25］。蒙特卡羅仿真時需要依靠大量光子才能得到較為準確的結果，在實際仿真時，考慮到計算機的內存空間使用情況，故只關注探測器處接收到的光子，對這些被接收到的光子路徑進行保存，如圖2所示。

Fig.2 Schematic diagram of light source and detector placement圖2 光源與探測器位置示意圖

光源和探測器放置的位置如圖2（a）所示，X 軸為光源和探測器所在的水平方向，Y 軸為頭皮表層的豎直方向，Z軸為腦組織深度方向。圖2（b）中光子入射點（0，0，0）與探測器的距離為SD，探測器與X 軸的夾角為θ。一定數量的光子從原點入射，在組織內進行隨機傳輸運動，部分光子可到達接收面并被接收。為便于后續處理，將探測器接收面看成邊長為Δs 的正方形。只有成功到達接收面的光子，才記錄其在組織內的傳輸路徑，其他光子則不予記錄。

仿真程序中所使用的五層顱腦模型光學參數如表1所示［26-28］。表中各參數的定義如下：n是折射率，g是各向異性因子，μa/cm-1是吸收系數，μs/cm-1是散射系數，t/cm 是各顱腦層的組織厚度。

Table 1 Optical parameters of the five-layer cranial model表1 五層顱腦模型光學參數

2 蒙特卡羅仿真與討論

為便于處理，本文將探測器與光子入射點設置為同一水平線，即設置圖2（b）中的θ值為0°，接收面的邊長Δs 為1mm，從原點（0，0，0）入射的光子數為107個，d的閾值為1mm。設置完成后，選擇不同的光源—探測器（SD）距離值運行MCML 程序。值得一提的是，由于此仿真需要對大量數據進行統計分析，故對電腦內存要求較高。仿真硬件條件為：Intel Core i9-12900處理器，128GB 的RAM 內存。

2.1 光子傳輸路徑

按如上條件設置后，可得到不同SD 距離下的光子傳輸路徑數據，選取SD=1cm 和3.5cm 的數據進行處理，繪制光子在組織內的傳輸路徑圖，如圖3所示。

Fig.3 Photon transport paths in tissue圖3 光子在組織內傳輸路徑

圖3（a）、（b）分別是SD=1cm 和SD=3.5cm 的傳輸路徑三維圖，中心區域是一團聚集在一起的密集路徑，邊緣區域的路徑逐漸稀疏。由圖3（a）和（b）可明顯看到，SD=3.5cm 時的傳輸路徑核心區域比SD=1cm 時大一些。

傳輸路徑圖反映了光子在顱腦組織內傳輸的部分情況，而確定光源—探測器（SD）間的最優距離，還需要借助其他參數進行判斷。本文通過研究光子傳輸時的最大可達深度和散射距離，作為選擇光源—探測器最優距離的參考。

2.2 光子最大可達深度

探測深度是腦血氧檢測中的一個重要參數，光子傳輸路徑的最大可達深度一定程度上反映了光子在顱腦內可探測的深度。根據每條光子傳輸路徑上Z 軸坐標的極大值，繪制出如圖4 所示的不同SD 距離下的光子最大可達深度概率分布圖。圖中橫坐標是光子傳輸的深度值，縱坐標是光子最大可達深度的概率值。可以看出，不同SD 距離對應的概率分布大致是正態分布。圖4（a）稍有不同，左側區域的概率值比右側區域密集些，主要是由于SD 的距離過小，光子多在淺層區域傳輸，傳輸到深層組織的光子相對較少。這也說明當SD<2cm 時，只有相對較少的光子可到達深層組織，這對于選擇近處探測器與光源之間的距離有很重要的參考價值。

Fig.4 Probability distribution of maximum reachable depth for 300 photonic transmission paths at different SD distances圖4 不同SD距離下300條光子傳輸路徑的最大可達深度概率分布

為便于分析，計算不同的SD 距離下光子最大可達深度在不同組織層的占比情況。如圖5 所示，橫坐標是不同的SD 距離，縱坐標是各顱腦組織在不同SD 距離時的最大深度占比情況。可看到隨著SD 距離的增大，頭皮層和顱骨層的占比不斷降低；腦脊液層的占比非常小，可忽略不計；灰質層的占比基本大于15%。頭皮和顱骨層占比在SD=1cm 處取得最大值23.78%。SD=2cm～4.5cm 時，光子最大可達深度在灰質層的占比于SD=3.5cm 處取得最大值19.59%

2.3 光子傳輸路徑在不同層的散射距離

僅采用最大可達深度描述光子對深層顱腦組織探測的情況還不夠客觀、全面，最大可達深度只能了解光子是否到達某一深度值，其他信息不得而知。

光子在顱腦組織中傳輸時，散射距離也是一個重要參數。因此，計算不同SD 距離下光子傳輸路徑在各層（取每一層深度為0.1cm）散射距離占路徑總長的比例，并繪制散射距離占比變化趨勢圖，如圖6所示。

Fig.6 Scattering distance variation of photon at different depths under different SD distances圖6 不同SD距離下光子在不同深度的散射距離變化

圖6 橫坐標表示顱腦組織的深度，縱坐標表示光子在各層的散射距離占傳輸路徑總長的比值。隨著深度增加，光子在深層組織處的散射距離所占比重不斷減少，表明光子在深層組織傳輸時的散射距離短。當深度達到3.5cm時，光子的散射距離所占比重基本為0，說明光子很難在大于3.5cm 的深層組織傳輸。通過圖6 可以了解光子散射距離的整體趨勢，但還不足以了解光子在不同顱腦層的情況。故繪制光子在不同顱腦層的散射距離占總散射距離的比重，如圖7所示。

Fig.7 The proportion of scattering distance of different brain layers to the tatal scattering distance圖7 不同顱腦層散射距離占總散射距離的比重

在圖7 中，橫坐標是不同的SD 距離，縱坐標是光子分別在頭皮和顱骨、腦脊液、灰質層的散射距離占總散射距離的比重。可明顯看到，光子的散射距離在頭皮和顱骨層所占比重很大，基本都在25%左右。由于腦脊液層的特殊光學性質（低吸收、低散射），故光子在腦脊液層的散射距離比重都很小。光子在灰質層的散射距離比重不等，與SD 距離相關，當SD=2～4.5cm 時，灰質層的散射距離比重在SD=4.5cm 處取最大值10.22%。這是因為隨著SD 距離的不斷增大，與較小的SD 相比，有更多光子可以到達更深層的組織進行傳輸。故隨著SD 的增大，灰質層的散射距離比重也在增大。

3 結果分析

根據多距離測量方法對SD 距離的選擇要求，近距離探測器D1 檢測頭皮和顱骨層的血液信息，遠距離探測器D2 檢測腦灰質層的血液信息。在對距離進行選擇時，需要同時考慮3 個方面：一是光源在傳輸時可達的最大深度；二是光源在組織中傳輸時攜帶信息多少；三是探測器實際接收到的信號強度。

有效深度比（Effective Depth Ratio，EDR）是探測器接收到來自有效深度層的光子與所有組織層光子的比值。將頭皮和顱骨層看作第一層，將腦脊液、灰質和白質層看作第二層，則在式（2）中，第一層（最大深度）為最大可達深度在頭皮和顱骨層，且并未穿過顱骨層的光子，第二層（最大深度）為穿透顱骨層的光子。

根據上文所述，計算出不同SD 距離下的EDR 值并繪制圖像，如圖8（a）所示。圖8（a）中，橫坐標為不同的SD 距離，縱坐標為EDR 值。從圖8（a）中可以看出，隨著SD 距離的增大，有效深度比EDR 逐漸減小，說明隨著SD 的增大，傳輸到深層組織的光子增多。EDR 值越大，說明SD 的距離越小，光子傳輸時位于頭皮和顱骨層的光子相對較多；反之，SD 距離越大，位于灰質層的光子相對較多。故可將EDR 值作為選擇光源與探測器間距離的重要依據。

Fig.8 Effective depth ratio and grey matter layer information ratio at different SD distances圖8 不同SD距離下的有效深度比與灰質層信息比

此外，光子在組織中傳輸時，攜帶的信息量大小與光子散射距離有很大關系。光子傳輸的EDR 值并不能完全反映探測器接收到來自灰質層的信息量大小，故對圖7 中的結果作進一步處理與分析。計算在不同SD 距離下光子傳輸時的灰質層信息比（Gray Matter Layer Information Ratio，GMLIR），表示為：

式中，灰質層（散射距離）為光子在灰質層的散射距離與總散射距離的比；頭皮、顱骨層（散射距離）為光子在頭皮和顱腦層的散射距離與總散射距離的比。繪制如圖8（b）所示的圖像，從圖中可明顯看出，隨著SD 距離的增大，灰質層信息不斷增多。在腦血氧檢測中，灰質層信息比的增大意味著探測器探測到的有效信息增多，這對于實際檢測是很重要的。故GMLIR 值也作為選擇光源與探測器距離的一個重要依據。

在實際選擇光源與探測器距離時，還需要考慮探測器接收到的信號強度。本文在研究中將光電探測器貼近頭皮表層進行必要的遮光處理后，分別使用760nm 和850nm的激光二極管（LD）在不同的SD 距離下照射頭皮表層，記錄探測器的電壓值。如圖9 所示，當不使用LD 照射時，探測器的電壓均值為3.67mV（參考值），而當使用LD 時，LD和探測器間的距離與信號強度成反比。SD 距離較小時，如SD=1cm，探測器接收到760nm 和850nm 信號的電壓均值大于400mV。隨著SD 距離的增大，信號強度開始明顯變小，說明SD 距離對于信號強度的影響非常大。探測器與LD 相距4cm 的位置處，探測器的電壓均值非常接近參考值，意味著探測器接收到的信號已很微弱。在確定SD距離時，需要考慮到探測器接收信號強弱這一方面的影響，故在選擇SD 時，將不再考慮SD>4cm 的情況。

Fig.9 Signal intensity received by the detector at different wavelengths with different SD distances圖9 不同SD距離時不同波長下探測器接收到的信號強度

根據以上分析，EDR 與GMLIR 是選擇SD 距離的重要依據。近距離探測器D1 主要檢測頭皮和顱骨層信息，需要EDR 值較大，且不檢測灰質層信息，因此需要GMLIR 值較小；遠距離探測器D2 主要檢測灰質層信息，其中含有頭皮和顱骨層的干擾信息少，故需要EDR 值較小，GMLIR 值較大。在選取光源與D2 的距離時，還需考慮探測器實際接收信號的強度。綜上所述，本文選擇的一組最優光源與探測器距離為SD1=1cm，SD2=3.5cm。

4 結語

本文通過選取一組顱腦組織的特性參數，利用蒙特卡羅方法對顱腦組織模型進行仿真研究，獲得光子在頭部組織的傳輸路徑，進而分析光子傳輸時的最大可達深度和散射距離，以及有效深度比EDR 和灰質層信息比GMLIR。根據多距離測量方法的要求，考慮了EDR、GMLIR 和探測器實際接收到的信號強度，確定一組適合選定顱腦組織參數的最優光源—探測器距離，SD1=1cm 和SD2=3.5cm。該組距離不僅可以最大程度地檢測灰質層的血氧信息，而且可以最大程度地處理來自頭皮和顱骨層的干擾，從而獲得純度更高的局部腦血氧信號。利用多距離測量方法設計的腦血氧檢測裝置，可以選擇此組最優距離進行傳感器探頭設計。本研究結果對腦血氧檢測裝置的設計，以及利用近紅外光譜技術檢測腦血氧的精度優化具有一定指導意義。