湯 偉,黃 璜,丑幸幸
(陜西科技大學(xué) 電氣與控制工程學(xué)院,陜西 西安 710021)
近年來(lái),人口老齡化發(fā)展形勢(shì)嚴(yán)峻,骨質(zhì)疏松癥[1]發(fā)病率逐年升高,現(xiàn)已成為人們生活中常見(jiàn)的代謝性疾病之一.人體骨骼在結(jié)構(gòu)上分為皮質(zhì)骨和松質(zhì)骨,而松質(zhì)骨由大量的骨小梁交錯(cuò)組織而成.骨質(zhì)疏松癥主要表現(xiàn)為松質(zhì)骨退化,骨密度下降,骨脆性增加[2].根據(jù)中國(guó)衛(wèi)生健康委員會(huì)的調(diào)查結(jié)果顯示,我國(guó)當(dāng)下存在骨密度檢測(cè)率低、低骨量人群年輕化、骨質(zhì)疏松癥認(rèn)知率低等問(wèn)題[3].臨床對(duì)于骨質(zhì)疏松癥的建議是預(yù)防大于治療,因此,及時(shí)的骨質(zhì)評(píng)價(jià)對(duì)預(yù)防骨質(zhì)疏松癥十分必要.
目前基于射線的骨質(zhì)檢測(cè)手段,如DXA、X線、定量CT 等輻射大,不適宜人群定期使用;而基于定量超聲的檢測(cè)方法無(wú)輻射、可靠性高,更適用于社區(qū)普查、孕婦復(fù)檢等情景.基于定量超聲對(duì)松質(zhì)骨及其骨質(zhì)疏松的評(píng)價(jià)方法主要分為透射法和背散射法[4],其中透射法將信號(hào)在骨中的衰減和聲速作為骨質(zhì)評(píng)價(jià)的主要參數(shù),但是很少反映骨微結(jié)構(gòu)信息,如骨小梁厚度、間距等,而這些微結(jié)構(gòu)信息已被證實(shí)與骨質(zhì)疏松有直接關(guān)聯(lián)[5].而背散射法通過(guò)骨小梁?jiǎn)卧裙俏⒔Y(jié)構(gòu)對(duì)超聲信號(hào)的散射作用,使接收的背散射信號(hào)包含了豐富的骨密度和骨微結(jié)構(gòu)信息.但是背散射法還需要建立更真實(shí)的骨模型和提取骨微結(jié)構(gòu)信息更準(zhǔn)確的算法[6].
在背散射信號(hào)分析算法發(fā)展中,多元線性回歸模型沒(méi)有考慮超聲參數(shù)與骨參數(shù)之間的非線性聯(lián)系[7],骨質(zhì)評(píng)價(jià)準(zhǔn)確度低;機(jī)器學(xué)習(xí)算法需要手動(dòng)提取特征信息,在處理大量數(shù)據(jù)時(shí)操作繁雜[8];而CNN 通過(guò)最小化預(yù)測(cè)誤差得到超聲參數(shù)與骨參數(shù)之間的非線性關(guān)系[9],可以自主進(jìn)行特征學(xué)習(xí),相比其他算法,評(píng)價(jià)結(jié)果更為準(zhǔn)確.
本文基于超聲背散射方法,首先對(duì)骨樣本進(jìn)行μ-CT 斷層掃描,獲得大量骨樣本原始圖像,對(duì)原始數(shù)據(jù)圖像進(jìn)行處理,導(dǎo)入Mimics中重建三維松質(zhì)骨模型;然后將重建好的不同結(jié)構(gòu)的骨模型導(dǎo)入Wave3000 Plus軟件中,模擬超聲在骨模型中的傳播過(guò)程,求解波動(dòng)方程,獲得大量背散射信號(hào);最后將獲得的背散射信號(hào)分為多組數(shù)據(jù)集輸入至改進(jìn)的CNN 網(wǎng)絡(luò)中訓(xùn)練測(cè)試,對(duì)骨樣本進(jìn)行預(yù)測(cè)評(píng)價(jià).
在骨骼結(jié)構(gòu)中,松質(zhì)骨雖然占骨量的20%,但構(gòu)成了80%的骨表面.現(xiàn)階段對(duì)松質(zhì)骨模型的研究,單圓柱模型[10]和多圓柱模型[11]與真實(shí)骨骼存在較大誤差,如表1所示.使用計(jì)算機(jī)仿真技術(shù)對(duì)骨科領(lǐng)域中的相關(guān)問(wèn)題進(jìn)行分析研究,避免了傳統(tǒng)離體和在體試驗(yàn)方法成本高、周期長(zhǎng)的問(wèn)題.Mimics能夠?qū)ΧS醫(yī)學(xué)影像信息進(jìn)行數(shù)字化三維重建,通過(guò)骨樣本CT 圖像數(shù)據(jù)三維重建出的松質(zhì)骨模型更近似真實(shí)骨骼.
表1 不同骨模型對(duì)比
本文對(duì)松質(zhì)骨的重建過(guò)程包括原始數(shù)據(jù)處理、數(shù)據(jù)導(dǎo)入、輪廓提取、區(qū)域增長(zhǎng)、模型建立.
(1)原始數(shù)據(jù)處理
本文原始圖像數(shù)據(jù)來(lái)自骨樣本斷層掃描,掃描儀型號(hào)Micro-CT Skyscan1076,儀器工作參數(shù)在電源75 KV,電流130μA,功率10 W,松質(zhì)骨樣本取自真實(shí)骨骼,在去除皮質(zhì)骨等預(yù)處理之后,根據(jù)骨的承重方向切割成表面平整的立方體結(jié)構(gòu)(約為20 mm*20 mm*20 mm).
對(duì)于原始數(shù)據(jù)圖像的預(yù)處理操作過(guò)程如下:首先μ-CT 對(duì)骨樣本掃描得到松質(zhì)骨CT 灰度圖像,建立松質(zhì)骨二值圖像.然后選取一幅松質(zhì)骨灰度圖像作為原始圖像,從任意方向隨機(jī)對(duì)骨小梁表面進(jìn)行腐蝕,將松質(zhì)骨二值模型中的l(骨小梁)替換為0(骨髓),模擬骨流失過(guò)程,松質(zhì)骨的孔隙度不斷增大.最后得到不同孔隙度結(jié)構(gòu)的松質(zhì)骨圖像模型,如圖1所示.
圖1 圖像預(yù)處理
(2)數(shù)據(jù)導(dǎo)入
Mimics支持DICOM、TIFF、BMP 等多種圖像格式,本文導(dǎo)入處理后的圖像數(shù)據(jù)為BMP 格式.進(jìn)入軟件界面后點(diǎn)擊File,打開(kāi)New Project,導(dǎo)入處理后的CT 斷層掃描圖像,層距為1 mm,在界面窗口中調(diào)整影像所示方向,其中冠狀位、矢狀位、橫斷位三個(gè)位置界面顯示的大寫(xiě)字母分別是:A、P和L、R 以及T、B,分別代表重建模型的前面、后面、左側(cè)、右側(cè)、上部、下部,可以通過(guò)上述不同視圖快速定位.
(3)輪廓提取
在操作窗口選擇閾值工具,設(shè)置骨組織CT值,閾值單位為HU,最小值是226 HU.在設(shè)置閾值時(shí),閾值既不能太低,也不能太高.閾值過(guò)低會(huì)受到軟組織干擾,過(guò)高則會(huì)造成骨組織信息丟失,本文設(shè)置閾值為255 HU 時(shí),分割效果最佳,如圖2所示.
圖2 輪廓分割
(4)區(qū)域增長(zhǎng)
選擇區(qū)域增長(zhǎng)按鈕,點(diǎn)擊感興趣區(qū)域,具有近似灰度值的骨骼周?chē)M織將被選中挑選除去,如圖3所示.
圖3 感興趣區(qū)域選取
(5)模型建立
選擇合適的蒙板,點(diǎn)擊Calculate重建按鈕,操作界面生成三維立體骨骼模型,可以通過(guò)平移、縮放、旋轉(zhuǎn),直觀化展現(xiàn)重建的松質(zhì)骨模型結(jié)構(gòu),如圖4所示.本文根據(jù)正常人群和骨質(zhì)疏松人群骨骼孔隙度的大小,重建了正常骨模型和異常骨質(zhì)疏松模型共四組,并導(dǎo)出為STL格式文件.
圖4 骨模型三維重建
Wave3000 Plus是一種超聲波仿真軟件,適用于無(wú)損測(cè)試材料評(píng)估和超聲醫(yī)學(xué)等研究領(lǐng)域.該軟件可以從STL文件中導(dǎo)入任何形狀的3D 對(duì)象,釆用時(shí)域差分方法[12](FDTD)求解完全的粘彈性波動(dòng)方程,求解快速,一次模型求解可以得到全場(chǎng)波的傳播,有利于大量數(shù)據(jù)的采集.FDTD通過(guò)對(duì)超聲背散射信號(hào)在松質(zhì)骨中的傳播進(jìn)行數(shù)值仿真,在考慮粘彈性材料對(duì)聲波吸收的基礎(chǔ)上,定義吸收邊界條件,模擬超聲在重建骨模型中的傳播過(guò)程.
其中,波動(dòng)方程形式如式(1)所示:
式(1)中:ρ代表材料密度,單位kg/m2;λ表示第一個(gè)Lame波常數(shù),單位N/m2;μ表示第二個(gè)Lame波常數(shù),單位N/m2;η表示剪切粘度,φ表示堆積粘度,單位N·s/m2;▽表示梯度算子;▽·表示散度算子;?表示偏微分算子;t表示時(shí)間,單位s;w可用式(2)計(jì)算得到.
仿真實(shí)驗(yàn)步驟為:首先導(dǎo)入STL 格式三維重建好的骨骼模型,3D 模型邊界按照空間坐標(biāo)分為Xmin、Xmax、Ymin、Ymax、Zmin、Zmax 六個(gè)平面;模型導(dǎo)入后,對(duì)模型材料進(jìn)行定義.
選取庫(kù)中材料骨骼和血液,骨模型樣本被假定浸在血液中,并假定血液填充所有孔隙空間,根據(jù)材料庫(kù)中材料固定參數(shù)如表2所示.
表2 骨模型材料參數(shù)
實(shí)驗(yàn)邊界條件設(shè)置為Xmin 平面為吸收層,Ymin、Ymax、Zmin、Zmax 平面為完全匹配吸收層,吸收層的作用是吸收傳播過(guò)程中的干擾信號(hào);然后定義發(fā)射接收源,背散射法主要通過(guò)一個(gè)超聲探頭對(duì)信號(hào)進(jìn)行發(fā)射接收,實(shí)驗(yàn)在Xmax平面設(shè)置一個(gè)1 MHz的正弦脈沖信號(hào)源,設(shè)置工作條件為Pluse echo方式.最后設(shè)置工作參數(shù)時(shí)間為30μs,步長(zhǎng)比例為0.8,分辨波長(zhǎng)為0.5 mm,周期為5,仿真結(jié)果如圖5所示.從仿真結(jié)果可以看出,波面被骨小梁嚴(yán)重散射.
圖5 背散射信號(hào)
圖6是一個(gè)背散射信號(hào)單周期時(shí)域波形示例.在信號(hào)起點(diǎn)之后,可以觀察到一個(gè)明顯的強(qiáng)反射回波,這是因?yàn)槌曉趥鞑ミ^(guò)程中遇到液體表面產(chǎn)生較大聲阻抗差異時(shí)出現(xiàn)的強(qiáng)反射,因此在獲得背散射信號(hào)后,需要對(duì)獲得的原始信號(hào)進(jìn)行處理.劉成成等[13,14]在選取背散射信號(hào)的有效區(qū)間時(shí),取了多組T1、T2的值,計(jì)算了不同背散射信號(hào)有效區(qū)間對(duì)應(yīng)的背散射參數(shù),通過(guò)分析背散射參數(shù)與骨密度之間的相關(guān)性,提出了背散射信號(hào)提取標(biāo)準(zhǔn).
圖6 背散射信號(hào)提取示例
該提取標(biāo)準(zhǔn)如下:從信號(hào)的起點(diǎn)位置去除時(shí)域部分時(shí)間長(zhǎng)度為T(mén)1的波段,以去除信號(hào)的末端做為新起點(diǎn),重新選取一段時(shí)間長(zhǎng)度為T(mén)2的信號(hào)作為有效信號(hào)區(qū)間(ROI),T2長(zhǎng)度包含的信號(hào)區(qū)間則是聲波在松質(zhì)骨模型中傳播的有效信號(hào).第五強(qiáng)強(qiáng)等[15]在該提取標(biāo)準(zhǔn)下取了相關(guān)性較高的背散射有效信號(hào)區(qū)間作為機(jī)器學(xué)習(xí)模型的輸入特征,骨質(zhì)評(píng)價(jià)實(shí)驗(yàn)的準(zhǔn)確率達(dá)到了82.86%.
本文在該提取標(biāo)準(zhǔn)下,對(duì)獲得的3 000組背散射信號(hào)進(jìn)行了相同的預(yù)處理.為了適應(yīng)后續(xù)CNN的輸入,需將預(yù)處理后的背散射信號(hào)數(shù)據(jù)進(jìn)一步進(jìn)行格式轉(zhuǎn)換,使背散射信號(hào)數(shù)據(jù)轉(zhuǎn)變?yōu)槎S數(shù)字矩陣,最后將轉(zhuǎn)換的數(shù)據(jù)切分成64*64尺寸大小的圖片,并制作成數(shù)據(jù)集.其次,將數(shù)據(jù)集中1 000組作為訓(xùn)練集,剩余2 000組劃分為四組測(cè)試集,訓(xùn)練集中包含有隨機(jī)數(shù)量正常骨模型[16]下仿真得到的背散射有效信號(hào).
(1)背散射信號(hào)分析方法
參考人群中骨質(zhì)疏松的臨床判斷標(biāo)準(zhǔn),通過(guò)背散射信號(hào)分析進(jìn)行骨質(zhì)評(píng)價(jià)時(shí),需要提取骨微結(jié)構(gòu)信息更準(zhǔn)確的算法.然而傳統(tǒng)的多元線性回歸模型并未考慮超聲參數(shù)與骨參數(shù)之間的非線性聯(lián)系,有監(jiān)督學(xué)習(xí)的超聲背散射方法和ANN[17]等機(jī)器學(xué)習(xí)算法需要手動(dòng)提取特征,操作繁雜.使用CNN可以將超聲參數(shù)映射到高維空間中,通過(guò)最小化預(yù)測(cè)誤差得到超聲參數(shù)與骨參數(shù)之間的非線性關(guān)系,在骨質(zhì)評(píng)價(jià)中優(yōu)于傳統(tǒng)算法.
ResNet作為經(jīng)典的深度CNN 模型,本身具有很好的非線性表達(dá)能力.但是經(jīng)過(guò)預(yù)處理和格式轉(zhuǎn)換后的背散射信號(hào)數(shù)據(jù)是二維的數(shù)字矩陣,只能從時(shí)間和空間維度提取一些特征信息來(lái)識(shí)別,會(huì)導(dǎo)致分類效果變差.因此針對(duì)使用數(shù)據(jù)的特性,對(duì)ResNet18模型進(jìn)行改進(jìn),提出一種ResNet16的改進(jìn)CNN 模型.劉銀萍等[18]已經(jīng)驗(yàn)證Res Net模型并不是所有的層都是必要的,較少的層結(jié)構(gòu)更有利于模型識(shí)別分類.在數(shù)據(jù)集數(shù)量不大的情況下,不需要多層次的殘差網(wǎng)絡(luò),否則易出現(xiàn)過(guò)擬合等現(xiàn)象[19].
(2)改進(jìn)網(wǎng)絡(luò)結(jié)構(gòu)
ResNet18的模型結(jié)構(gòu)包含了1 個(gè)卷積層,4個(gè)殘差塊結(jié)構(gòu),每個(gè)殘差塊包括4個(gè)卷積層,以及1個(gè)全連接層,通過(guò)平均池化與全連接層相連,池化的結(jié)果會(huì)使得特征減少和模型參數(shù)減少[20].平均池化是對(duì)鄰域內(nèi)特征點(diǎn)進(jìn)行平均運(yùn)算,而最大池化是對(duì)鄰域內(nèi)特征點(diǎn)取最大值.在CNN 模型中,特征提取的誤差主要來(lái)自兩個(gè)方面:(1)鄰域大小受限造成的估計(jì)值方差增大;(2)卷積層參數(shù)誤差造成估計(jì)均值的偏移.因此本文在擴(kuò)寬網(wǎng)絡(luò)結(jié)構(gòu)時(shí),使用平均池化和最大池化分別將深度特征信息進(jìn)行融合,有利于提高模型原始數(shù)據(jù)的特征信息利用率.若去除其余卷積層,雖然模型的參數(shù)會(huì)得到減少,但是特征信息卻不能得到有效融合,會(huì)使特征提取產(chǎn)生誤差[21].
本文改進(jìn)的網(wǎng)絡(luò)結(jié)構(gòu)參數(shù)與未改進(jìn)的Res-Net18和具有相同層數(shù)的VGG16結(jié)構(gòu)參數(shù)如表3所示.其中卷積層conv1的主要作用是通過(guò)卷積運(yùn)算減少輸入數(shù)據(jù)大小,殘差塊conv2_x、conv3_x、conv4_x、conv5_x 由多個(gè)卷積層組成,其作用是對(duì)縮減后的數(shù)據(jù)進(jìn)行特征提取.本文首先用背散射信號(hào)數(shù)據(jù)集從頭開(kāi)始訓(xùn)練;更改輸入數(shù)據(jù)的大小,減少訓(xùn)練時(shí)間;縮減ResNet18 網(wǎng)絡(luò)的深度,降低冗余參數(shù),去除一個(gè)殘差結(jié)構(gòu)conv5_x,分解為兩個(gè)并列的卷積層,拓寬整體模型的寬度,所用卷積核都是3*3大小,將conv5_x卷積層分成兩個(gè)并列的卷積層有兩個(gè)好處:首先拓寬CNN 模型的寬度有利于在提取深度特征信息時(shí)提高信息的利用率,進(jìn)而提高模型的準(zhǔn)確度;其次分解過(guò)后,使冗余參數(shù)得到了減少,降低了訓(xùn)練時(shí)間.分解之后的卷積層用最大池化和平均池化可以將得到的特征信息進(jìn)行深度融合,使CNN 模型提取的特征信息更加準(zhǔn)確[22].
表3 CNN 模型參數(shù)
改進(jìn)CNN 模型結(jié)構(gòu)如圖7 所示.由于CNN模型的訓(xùn)練對(duì)實(shí)驗(yàn)配置要求較高,本文實(shí)驗(yàn)的訓(xùn)練環(huán)境采用的是Windows10系統(tǒng),使用顯卡的型號(hào)為RTX2060,其中計(jì)算機(jī)最大內(nèi)存容量為16 GB,計(jì)算機(jī)CPU 型號(hào)為AMD R7-4800H,為八核十六線程處理器.模型的訓(xùn)練和編程軟件使用Python3.7版本和Anaconda3,環(huán)境后端框架為T(mén)ensor Flow.
圖7 改進(jìn)網(wǎng)絡(luò)結(jié)構(gòu)
最后將實(shí)現(xiàn)的16層ResNet16在數(shù)據(jù)集上訓(xùn)練測(cè)試,并與未改進(jìn)的Res Net18以及具有相同層數(shù)的VGG16進(jìn)行簡(jiǎn)單的性能比較.如圖8 所示,改進(jìn)模型訓(xùn)練效果的誤差率對(duì)比未改進(jìn)的Res-Net18和VGG16有明顯降低.
圖8 不同模型訓(xùn)練對(duì)比
由于本實(shí)驗(yàn)是為了對(duì)骨質(zhì)狀況進(jìn)行分類,所以使用的評(píng)價(jià)指標(biāo)主要為準(zhǔn)確率(Accuracy)和特異性(Specificity),其中準(zhǔn)確率可以由公式(3)和公式(4)來(lái)表示:
其中,TP(True Positive)代表真正例,即被正確預(yù)測(cè)為正例的個(gè)數(shù);FN代表假負(fù)例,即被錯(cuò)誤預(yù)測(cè)為負(fù)例的個(gè)數(shù);FP代表假正例,即被錯(cuò)預(yù)測(cè)為負(fù)例的個(gè)數(shù);TN代表真負(fù)例,即被正確預(yù)測(cè)為負(fù)例的個(gè)數(shù).
本文使用改進(jìn)的CNN 模型對(duì)測(cè)試集進(jìn)行測(cè)試,在迭代500次訓(xùn)練后,改進(jìn)CNN 模型收斂效果趨于穩(wěn)定.
將重建的不同骨模型下獲得的四組數(shù)據(jù)集作為輸入,對(duì)CNN 模型進(jìn)行測(cè)試,實(shí)驗(yàn)準(zhǔn)確率達(dá)到了97.3%,其中損失函數(shù)變化和準(zhǔn)確率變化如圖9所示.本文改進(jìn)的ResNet16模型的預(yù)測(cè)分類準(zhǔn)確率對(duì)比文獻(xiàn)15 中有監(jiān)督學(xué)習(xí)算法和文獻(xiàn)17 中ANN 算法分類準(zhǔn)確率結(jié)果如表4所示.
圖9 改進(jìn)模型準(zhǔn)確率和損失函數(shù)
表4 不同算法性能對(duì)比
通過(guò)不同算法之間的對(duì)比,實(shí)驗(yàn)結(jié)果表明ResNet18中的層可以進(jìn)行改進(jìn),而相同層數(shù)結(jié)構(gòu)下,ResNet模型基于跳層連接思想比VGG16 網(wǎng)絡(luò)和ANN 等算法能更好的表征特征信息,有利于深層特征的提取,提升模型分類準(zhǔn)確率.本文提出的改進(jìn)CNN 算法在提取骨微結(jié)構(gòu)信息上魯棒性更好,準(zhǔn)確度更高,對(duì)骨樣本有更好的區(qū)分能力.
本文在重建了不同孔隙結(jié)構(gòu)的松質(zhì)骨模型后,對(duì)超聲在骨模型中的背散射過(guò)程進(jìn)行了仿真,獲得了大量背散射信號(hào),解決了數(shù)據(jù)來(lái)源少的問(wèn)題.在分析背散射信號(hào)進(jìn)行骨質(zhì)評(píng)價(jià)時(shí)提出了一種基于深度學(xué)習(xí)的超聲背散射方法使用改進(jìn)的CNN 模型算法對(duì)骨樣本進(jìn)行預(yù)測(cè)分類.
研究結(jié)果表明:基于深度學(xué)習(xí)的超聲背散射方法對(duì)骨樣本的分類準(zhǔn)確率達(dá)到了97.3%,提取骨微結(jié)構(gòu)信息更加準(zhǔn)確,對(duì)骨質(zhì)疏松樣本有較強(qiáng)的區(qū)分能力.基于深度學(xué)習(xí)的超聲背散射評(píng)價(jià)效果優(yōu)于現(xiàn)有的其它定量超聲方法,骨質(zhì)評(píng)價(jià)更為可靠,對(duì)超聲背散射在體應(yīng)用具有參考意義.
本文的研究結(jié)果對(duì)超聲背散射在臨床應(yīng)用、智慧醫(yī)療、輔助診斷等方面有一定幫助,在后續(xù)模型優(yōu)化中,可以進(jìn)一步考慮結(jié)合群體信息提升模型的準(zhǔn)確率.本次實(shí)驗(yàn)由于無(wú)法得到離體樣本的群體信息,未能研究群體信息對(duì)骨質(zhì)疏松癥預(yù)測(cè)的貢獻(xiàn).