基于錐形束CT數據估算頜骨骨量變化的方法研究*

李偉志 王瑞永 劉 敏 顧曉宇 喬惠婷

頜骨的骨量變化，尤其是牙槽突骨量減少，是口腔醫學臨床工作中常常要面對的問題。牙周病的進展和治療會使牙槽突產生不同程度的骨量變化[1]。牙拔除術不可避免地導致牙槽骨骨量減少[2]，即使在牙齒拔除（或脫落）后自然愈合期牙槽突骨量減少還會進行性發展[3，4]。有效定量定位地評估頜骨骨量變化程度對于評價牙拔除術對牙槽骨的損傷程度、牙周病進程及正畸正頜美容整形等手術效果有重要意義[5]。近年來微創拔牙相關手術技術不斷發展，旨在減少對牙槽骨的損傷[6]；組織工程技術的應用（比如牙種植或牙周相關植骨手術）在一定程度上可以修復牙槽突骨量減少的情況[7]，牙槽突骨量變化的準確定位和對應體積的有效估算對于口腔臨床工作中新技術及新方法有重要的支撐意義[8]。

借助錐形束CT（cone-beam computed tomography，CBCT）可以無創地進行頜骨（含牙槽突）的骨量探測[9]。目前對于頜骨牙槽突骨量變化的檢測，主要通過人工設定標志點并進行測量完成[10，11]，這種方法對成像方位有嚴格要求[12]，得到的測量值僅為部分二維平面特征，無法對三維空間內頜骨牙槽骨的骨量變化進行全面評估，并且所獲得頜骨牙槽骨特征值受閱片人為因素影響較大。Linderup 等人[13]提出基于兩次CBCT 成像，分別建立頜面部三維模型，通過比較三維模型體積的變化來進行頜骨牙槽骨-骨量變化體積估算。由于涉及頜面部的三維建模工作量較大，不易在臨床使用，因此臨床急需自動或半自動的頜骨牙槽骨骨量變化的估算方法，來實現不同時期頜骨牙槽骨骨量變化情況的定量評估。

本文提出了一種基于兩次CBCT 成像數據進行半自動頜骨牙槽骨骨量變化定量定位估算的方法，通過在不同時間點（如牙拔除術前后）分別進行CBCT 成像，以初始狀態時刻的CBCT 圖像作為基準，利用最大化互信息配準法，人工結合計算機程序化算法半自動將骨量變化估計時刻的CBCT 圖像與基準進行配準，對比兩幅已配準的CBCT 圖像中頜骨差異，分離存在差異的體素，并根據成像參數估算兩次成像之間的頜骨牙槽骨的骨量變化。本文利用平臺實驗驗證該方法的可行性，并結合光學掃描實驗對該方法進行了量化校正，證明該方法可以明確定位頜骨骨量變化位置并實現頜骨骨量變化體積的準確估算。

1.方法描述

1.1 頜骨牙槽骨骨量（體積）變化估算方法流程

估算從t1 到t2 一段時間內頜骨牙槽骨骨量變化的方法，需要實施2 次CBCT 成像，在配準兩幅三維圖像后，自動實現圖像的差異提取及體素估算。首先，要在初始狀態t1 時刻對頜面部進行CBCT 成像，得到參考圖像，而后在擬進行頜骨牙槽骨骨量變化估算的t2 時刻再次對頜面部進行CBCT 成像，得到浮動圖像?；诠趋赖膭傂蕴攸c，配準算法選用剛體變換。估算剛體變換參數，把變換參數的初始值設為零。配準過程中利用估算參數對浮動圖像進行變換。選用三線性插值方法進行體素插值，使參考圖像與浮動圖像具有相同的空間分辨率。定義互信息為兩圖像的相似性測度，用它衡量兩幅圖像配準后的匹配程度。采用方向加速法，迭代執行上述步驟，通過調整幾何變換參數，使浮動圖像和參考圖像的互信息測度最大。配準完成后，參考圖像和浮動圖像做差獲得差異圖像。針對差異圖像進行閾值分割獲得頜骨牙槽骨變化量，完成一段時間內骨量實際變化的定量估算。

圖1 估算頜骨牙槽骨骨量變化的方法流程圖

1.2 基于最大互信息的圖像配準 參考圖像與浮動圖像的配準是自動或半自動實現缺損體積估計的重要基礎，可以克服兩次獨立成像體位不同等成像因素的影響，實現影像數據的有效對比，具體算法流程如圖2 所示?；バ畔⑹且粋€統計相關性的測度，通常用于表示兩個隨機變量的依賴程度，也可以表示兩幅圖像的匹配程度，兩幅圖像的互信息值最大時，配準精度最高[14]。參考圖像A 與浮動圖像B的互信息MI(A，B)定義如下：

其中，PA(a)和PB(b)分別是參考圖像A 和浮動圖像B 的灰度邊緣概率分布，PAB(a，b)為參考圖像A 和浮動圖像B 的聯合概率分布。直接計算兩幅圖像的聯合概率分布比較困難，可以借助于兩幅圖像的聯合直方圖h(a，b)獲得[15]。

圖像的精細配準過程是基于圖像互信息的目標函數優化問題，其目標是找到使浮動圖與參考圖像對齊的空間映射，即變換矩陣T（包含旋轉信息與位移信息）。變換矩陣可以實現圖像在三維空間的方位變換，本算法采用Powell 搜索求解最優化變換矩陣。三維空間優化過程是由4 次一維搜索組成，分別沿x，y，z 三個正交方向分別搜索得到單向最優值x’，y’和z’，再沿初始點與單向最優值對應點（x’，y’，z’）連線方法進行一維搜索得到最優值，多次迭代后可得到最優解。

圖2 基于最大互信息的配準算法流程圖

1.3 差異圖像的分割 將配準后的浮動圖像和參考圖像的對應體素作差后，獲得差異圖像。差異圖像可以去除牙槽骨周圍組織的影響，突出牙槽骨骨量變化。針對差異圖像再進行閾值分割，獲得兩次成像牙槽骨變化的區域。閾值分割是一種直接檢測區域的分割方法，它假設目標內相鄰像素間的灰度值是相似的，但目標與背景的像素在灰度值上有差異。在灰度直方圖上，目標和背景對應不同的峰值，選取的閾值位于兩個峰谷間，從而將目標和背景分開。由于受不同成像因素的影響兩次成像配準后并非完全相同，做差后得到的差異圖像會體現出成像差異，不同于組織無變化區域的噪聲數據，牙槽骨骨量變化區域會呈現明顯數值，因此對于配準后作差得到的差異圖像進行閾值分割可以選出三維圖像中骨量變化的體素。

1.4 頜骨牙槽骨變化量估算 在實現骨量變化區域進行精確分割后，針對已標記的變化體素，使用如下公式進行計算可以得到骨量變化的實際體積。頜骨牙槽骨的骨量變化體積的計算公式為：

其中：V 為估算出的骨量變化體積，N 表示該區域的斷層數目；Ai表示第i 斷層內-發生骨量變化的區域的面積；Δh表示CBCT 重建數據的斷層層厚，Mi為對差異圖像分割后第i 層所分割出的像素數，a0為每一斷層內單位像素所對應的面積，M 為對差異圖像分割后得到骨量變化區域的體素數量，V0為三維影像中單位體素所對應的面積。該實驗中所用CBCT 重建層厚0.25 mm，每一斷層內像素對應面積0.0625 mm2，每個體素體積為0.0156 mm3。

2.方法的實驗驗證

2.1 實驗對象 本文利用人類濕下頜骨標本及密度與人體牙槽骨接近的豬肋骨塊進行實驗，其中涉及1套人類濕下頜骨標本，該下頜標本存在10顆完整牙齒，5 顆殘根及殘冠，如圖3（1）所示，豬骨塊共8 塊，取自成年豬第4 肋，骨塊由小到大排列，如圖3（2）所示。豬肋骨塊與人體骨松質具有相近的成像特性，而發生在骨松質的下頜骨骨量改變非常不容易被察覺、辨別，因此本文通過在下頜骨標本右下頜第二磨牙遠中的磨牙后墊所對應的骨組織表面放置豬肋骨塊模擬多種下頜骨骨量變化情況。本研究獲得北京市海淀醫院生物倫理委員會審查批準（批準號2019022）。

圖3 CBCT成像實驗所用材料：（1）人體濕下頜骨，（2）豬骨塊

2.2 實驗過程 使用口腔錐形束CT（北京朗視儀器有限公司），對人類濕下頜骨標本及豬肋骨塊進行成像，其中對人類濕下頜骨標本單獨進行成像，并將8塊豬骨塊分別放置在人類濕下頜骨標本上，分別進行全方位掃描成像。所有研究對象掃描參數設定為：電流4 mA，電壓100 kV，將掃描數據用錐形束CT自帶的軟件進行三維圖像重建，重建后圖像以DICOM格式存儲，共得到三維錐形束CT圖像9套，其中1套為獨立人類濕下頜骨三維影像，另8套分別為帶有不同體積的豬骨塊的人類濕下頜骨三維影像，各體素均為0.25 mm×0.25 mm×0.25 mm。

獲得錐形束CT三維圖像后，基于VS2017平臺C++語言，在主頻2.2 GHz，Intel-Xeon E5-2650 v4處理器，64GByte內存環境下運行算法程序，分別對8組對照圖像進行配準、做差、分割后實現下頜骨骨量變化體積的估算。

鑒于CBCT存在容積效應，基于CBCT圖像通過算法的頜骨骨量缺損的估算值可能會出現少量偏差。本文使用光學三維掃描儀（Breuckmann precision in 3D），逐一對豬肋骨進行光學掃描得到豬肋骨塊的實際體積。并利用MATLAB 軟件，以光學掃描得到的體積為真值，通過一元線性回歸的方法對基于CBCT 的頜骨骨量缺損算法進行定量線性校正，為確保在頜骨骨量無變化時能準確估算，回歸方程常數項設為0，回歸系數即為算法的校正系數。

3.結果

本文分別對人體濕下頜骨標本以及放置體積不同的豬肋骨塊后的下頜骨標本進行CBCT 成像，共得到9套CBCT 圖像，形成8組對照圖像（放置豬骨塊的下頜骨標本與無豬骨塊下頜骨標本的2 次對照成像），圖4 是在使用CBCT 常用設置，未進行人工處理條件下得到CBCT 斷層圖像，這是一組對照圖像中的某一斷層，其中紅色圈出豬肋骨塊放置的位置。

圖4 人體濕下頜骨及豬骨塊CBCT成像

表1 基于圖像的骨量變化估算值與光學掃描值（mm3）

成像后，分別對8 組對照圖像進行配準、做差、分割后實現下頜骨骨量變化體積的估算，得到的骨量變化數據（見表1）。算法能將人眼不易察覺的骨量變化情況清晰地顯示出來，如圖5 所示，圖5（1）、（2）和（3）是利用CBCT（朗視）采用Ray-casting方法進行體繪制所得到的三維影像，其中圖5（1）、（2）均是下頜骨上放置豬肋骨塊后所得到的三維影像，在口腔科CBCT 常規設置下直接得到的三維圖像圖5（1）中并不能清晰地觀測豬肋骨塊，在人工手動調整顯示窗位、窗寬后圖5（2）中隱約可見豬肋骨塊，圖5（3）是無豬肋骨塊的下頜骨獨立成像所得到的對照三維影像，而圖5（4）是利用算法估算出頜骨骨量變化情況的可視化顯示，紅色區域為骨量變化區域，與成像過程中所放置的豬骨塊位置一致。

圖5 CBCT圖像三維可視化結果及算法標注骨量變化區域

通過光學掃描分別獲得8 個大小不同豬骨塊的實際體積（表1）作為骨量變化的真值對算法進行校正。所得到的校正系數為0.8044，校正后的確定系數R-square為0.988，校正結果如圖6所示。

圖6 骨量變化估算的定量校正結果

4.討論

在臨床中頜骨的骨量變化是難以直接測量的，因此頜骨的缺損與修復都難以進行定位定量分析。將同一對象的兩次錐形束CT成像與圖像處理算法相結合，應用于頜骨骨量變化的估算給口腔手術及治療技術的評價帶來了新的思路。本文提出基于兩次錐形束CT 成像數據估算頜骨骨量的變化值，并形成算法程序，可以通過牙拔除術前后兩次成像判斷牙拔除術對下頜牙槽骨的即刻損傷，也可以對牙周、正畸或正頜手術等治療后患者進行定期成像，跟蹤治療進展及此過程中頜骨骨量的長期變化，對于口腔科治療技術效果的定量評價有重要意義，可以證明并推動新治療技術的發展。

本文基于兩次錐形束CT成像，提出了半自動化估算頜骨變化的方法的具體流程，其中兩次獨立成像的錐形束CT圖像的配準是重要基礎，本文利用最大互信息原理實現了圖像的自動配準算法，并基于配準結果實現了做差、分割及骨量變化的定量定位估算算法。盡管錐形束CT成像受體位、具體操作設置及多種因素影響，任意兩次成像都存在差異，但基于互信息的配準算法有效避免了成像體位差異的影響因素。配準后做差得到的差異圖像大大降低了周圍組織的影響，將骨量變化凸顯出來。后續在差異圖像的感興趣區域進行閾值分割，通過閾值篩選去掉無變化區域的影響，將骨量變化的區域分割出來，通過在計算機上運行程序，導入兩次成像得到的CBCT圖像數據，可同時實現骨量變化的定位標記與定量估計。

為驗證方法的有效性，本文設計了利用人體濕下頜骨及豬肋骨塊的下頜骨改變仿真成像實驗，其中利用豬骨塊模擬下頜骨松質骨骨量變化。本文選取了下頜牙槽骨骨松質區域，與豬肋骨塊進行區域體素CT灰度值比較，人工選取的牙槽骨骨松質區域與豬肋骨塊的CT灰度值非常接近，平均相對差異為1.1%。在此方法中沒有其他人體組織CT灰度值介于豬骨和人骨之間對計算產生干擾，因此利用豬肋骨塊進行下頜骨骨量變化成像實驗是有效并且可行的。

基于8 組成像實驗，對所提出方法進行了校驗。對于多次獨立成像所得到了錐形束CT頜骨圖像，本文所提出的基于最大互信息的配準方法可以在不做特殊標記的前提下自動實現圖像的配準，并估算出頜骨骨量變化并顯示出相關變化的空間對應位置。隨著成像所用豬肋骨大小的變化，算法能很好地估算出下頜骨骨量變化情況。利用光學掃描豬肋骨塊得到每組成像對象的真實骨量差異，對該方法的估算結果進行比對，發現一致性很好，隨著真實骨量差異增大，算法估算結果隨之增大，并且算法估算結果與真實值具有很好的線性關系。這主要是由于CT 成像存在容積效應，基于CT 圖像所估算的體積值會與真實體積有偏差。本文利用光學掃描結果對算法估算值進行修正，得到修正系數0.8044，修正后確定系數為0.988。確定系數r-square 的正常取值范圍為[0，1]，其越接近1 表面估算結果與真值越接近，因此校正后的算法可以對0-500 mm3間骨量變化進行有效估算。

本文提出的方法是基于錐形束CT圖像進行骨量變化的估算，估算方法的靈敏度受限于成像體素大小及算法重采樣體素大小，該方法對于骨量變化估算的靈敏度為0.0156 mm3，高于人眼的靈敏度。該方法不僅可以定量估算出頜骨骨量變化，還可以標記出骨量變化的位置。如圖4中紅色顯示區域，在算法中具有做差分割功能，可以將骨量變化區域單獨分割出來，若與原圖像融合即可清晰地顯示頜骨骨量變化的區域，這對于口腔科手術及治療的相關研究有重要價值。定量的骨量變化估算是評價手術等治療手段效果的基本指標之一，而同時明確骨量變化的位置對于技術手段的改進有積極的指導意義。