基于SE-DenseVoxNet的骨小梁模量預(yù)測(cè)方法

曹熠煒，何思淵，周平，李瀾，包倪榮

1 東南大學(xué) 生物科學(xué)與醫(yī)學(xué)工程學(xué)院，南京市，210096

2 南京大學(xué)醫(yī)學(xué)院附屬鼓樓醫(yī)院，運(yùn)動(dòng)醫(yī)學(xué)與成人重建外科，醫(yī)藥生物技術(shù)國(guó)家重點(diǎn)實(shí)驗(yàn)室，南京市，210008

3 解放軍東部戰(zhàn)區(qū)總醫(yī)院，南京市，210002

0 引言

骨質(zhì)疏松癥是一種常見(jiàn)的全身代謝性疾病，具有致殘率高、易發(fā)骨折等特點(diǎn)。隨著人口老齡化進(jìn)程的加速，骨質(zhì)疏松癥的發(fā)病率逐漸上升，嚴(yán)重影響老年人的健康和生活質(zhì)量[1]。骨質(zhì)疏松癥患者體內(nèi)的骨吸收和骨形成的動(dòng)態(tài)平衡被破壞，使得機(jī)體的骨質(zhì)流失，骨小梁結(jié)構(gòu)退化，骨的生物力學(xué)性質(zhì)下降，容易引發(fā)骨質(zhì)疏松性骨折[2]。因此，準(zhǔn)確預(yù)測(cè)骨的生物力學(xué)性質(zhì)對(duì)于骨質(zhì)疏松性癥的早期防治與診斷具有重要的意義。

目前國(guó)內(nèi)外主要通過(guò)雙能X射線(xiàn)法(DXA)測(cè)量骨骼的骨密度(BMD)[3]，進(jìn)而評(píng)估和預(yù)測(cè)骨的生物力學(xué)性質(zhì)，對(duì)骨質(zhì)疏松癥進(jìn)行診斷。這種方法存在兩個(gè)問(wèn)題：第一，雙能X射線(xiàn)法測(cè)量骨密度時(shí)會(huì)受到皮質(zhì)骨外殼、脂肪組織和軟組織的干擾，測(cè)量結(jié)果存在誤差；第二，骨密度是一個(gè)宏觀參數(shù)，無(wú)法完全反映骨骼的形態(tài)結(jié)構(gòu)特征。研究表明，骨密度只能反映60%～70%的骨強(qiáng)度變化，骨的生物力學(xué)性質(zhì)與骨小梁微結(jié)構(gòu)密切相關(guān)[4]。Micro-CT、MRI等三維醫(yī)學(xué)成像方法可以消除骨骼外圍組織的干擾，獲取骨小梁微結(jié)構(gòu)的細(xì)節(jié)，被廣泛應(yīng)用于骨骼樣本的觀測(cè)與研究。ULRICH等[5]基于Micro-CT圖像對(duì)骨骼樣本的骨體積分?jǐn)?shù)、各向異性程度、結(jié)構(gòu)模型指數(shù)、骨小梁間隙等參數(shù)進(jìn)行統(tǒng)計(jì)，然后使用多元線(xiàn)性回歸的方法對(duì)樣本的彈性模量進(jìn)行預(yù)測(cè)。實(shí)驗(yàn)結(jié)果表明，多元線(xiàn)性回歸方法的預(yù)測(cè)精度優(yōu)于骨密度的預(yù)測(cè)精度；XU等[6]和HUBER等[7]分別基于Micro-CT圖像和MRI圖像統(tǒng)計(jì)了樣本的骨體積分?jǐn)?shù)、面積體積比、骨小梁數(shù)目等參數(shù)，然后使用SVR、kNN等機(jī)器學(xué)習(xí)算法對(duì)樣本的屈服應(yīng)力進(jìn)行預(yù)測(cè)。實(shí)驗(yàn)證明，機(jī)器學(xué)習(xí)算法的預(yù)測(cè)誤差小于多元線(xiàn)性回歸算法；CIUSDEL等[8]通過(guò)卷積神經(jīng)網(wǎng)絡(luò)的方法自動(dòng)提取骨小梁的形態(tài)結(jié)構(gòu)特征，對(duì)樣本在指定載荷下的應(yīng)變進(jìn)行預(yù)測(cè)。實(shí)驗(yàn)結(jié)果表明，卷積神經(jīng)網(wǎng)絡(luò)的預(yù)測(cè)精度相比機(jī)器學(xué)習(xí)方法有所提升。

為了更準(zhǔn)確地預(yù)測(cè)骨的生物力學(xué)性質(zhì)，嘗試使用深度學(xué)習(xí)方法自動(dòng)提取骨小梁微結(jié)構(gòu)的形態(tài)特征，對(duì)樣本的彈性模量進(jìn)行預(yù)測(cè)。對(duì)DenseVoxNet[9]的網(wǎng)絡(luò)結(jié)構(gòu)進(jìn)行改進(jìn)，提出了深度回歸網(wǎng)絡(luò)SE-DenseVoxNet，并基于骨小梁三維圖像數(shù)據(jù)集對(duì)模型進(jìn)行了訓(xùn)練和評(píng)估。

1 方法

DenseVoxNet是一種三維圖像分割網(wǎng)絡(luò)，由特征提取和上采樣兩部分組成。我們對(duì)DenseVoxNet的特征提取部分進(jìn)行改進(jìn)，提出了深度回歸網(wǎng)絡(luò)SE-DenseVoxNet。如圖1所示，SE-DenseVoxNet首先基于卷積層Conv1對(duì)輸入圖像降維，然后通過(guò)Dense模塊[10]和SE模塊[11]提取圖像特征，最終通過(guò)全局平均池化層和全連接層對(duì)樣本真值進(jìn)行預(yù)測(cè)。如圖2所示，Dense模塊由12層網(wǎng)絡(luò)層密集連接構(gòu)成，每層網(wǎng)絡(luò)層的結(jié)構(gòu)為歸一化層-ReLU函數(shù)-卷積層-歸一化層-ReLU函數(shù)-卷積層，其中第一層卷積層包含48個(gè)1×1×1大小的卷積核，第二層卷積層包含12個(gè)3×3×3大小的卷積核。1×1×1的卷積操作降低了特征維度，有效地減少了模型的參數(shù)數(shù)目和計(jì)算量。通過(guò)網(wǎng)絡(luò)層之間的密集連接，Dense模塊實(shí)現(xiàn)了特征復(fù)用，改善了網(wǎng)絡(luò)的信息流和梯度流，使得深層網(wǎng)絡(luò)的訓(xùn)練變得更容易。

圖1 SE-DenseVoxNet網(wǎng)絡(luò)結(jié)構(gòu)示意圖Fig.1 The illustration of proposed DenseVoxNet

圖2 Dense模塊結(jié)構(gòu)示意圖Fig.2 The illustration of Dense module

SE模塊的網(wǎng)絡(luò)結(jié)構(gòu)，如圖3所示。SE模塊分為左右兩路結(jié)構(gòu)，左路結(jié)構(gòu)首先通過(guò)1×1×1的卷積操作對(duì)輸入特征的通道維度進(jìn)行壓縮，生成各個(gè)空間位置的權(quán)重，然后基于空間位置的權(quán)重對(duì)輸入特征進(jìn)行標(biāo)定；右路結(jié)構(gòu)首先通過(guò)全局平均池化層對(duì)輸入特征的空間維度進(jìn)行壓縮，然后通過(guò)兩層全連接層生成了各個(gè)通道的權(quán)重，最終基于各個(gè)通道的權(quán)重對(duì)輸入特征進(jìn)行標(biāo)定。SE模塊匯總了左右兩路結(jié)構(gòu)的輸出，從空間維度和通道維度上對(duì)輸入特征進(jìn)行了自適應(yīng)的校準(zhǔn)，增強(qiáng)了有用的特征，抑制無(wú)效的特征。

在第一個(gè)SE模塊和第二個(gè)Dense模塊之間插入了過(guò)渡層(transition layer)[10]，以減少模型的參數(shù)數(shù)目。過(guò)渡層由歸一化層、ReLU函數(shù)、卷積層Conv14和最大池化層構(gòu)成，其中卷積層Conv14包含160個(gè)1×1×1大小的卷積核，用于整合特征通道之間的信息。最大池化層的步長(zhǎng)為2，將特征維度降低為原來(lái)的一半。SEDenseVoxNet的整體網(wǎng)絡(luò)結(jié)構(gòu)及參數(shù)如圖3和表1所示。

圖3 SE模塊結(jié)構(gòu)示意圖Fig.3 The illustration of SE module

表1 SE-DenseVoxNet的網(wǎng)絡(luò)結(jié)構(gòu)及參數(shù)Tab.1 The architecture and parameters of SE-DenseVoxNet

2 數(shù)據(jù)集

與南京市鼓樓醫(yī)院合作，收集了9個(gè)股骨頭壞死患者全髖關(guān)節(jié)置換術(shù)時(shí)取出的股骨頭樣本，然后使用Micro-CT對(duì)樣本進(jìn)行掃描重建，重建分辨率為78 μm，各向同性。通過(guò)感興趣區(qū)域提取、閾值分割、圖像切割、數(shù)據(jù)增強(qiáng)等步驟對(duì)樣本的Micro-CT圖像進(jìn)行預(yù)處理，共計(jì)得到249 864個(gè)32×32×32大小的骨小梁三維二值圖像，作為骨小梁三維圖像數(shù)據(jù)集的樣本。

得到了骨小梁三維圖像數(shù)據(jù)集的樣本后，使用有限元軟件Abaqus對(duì)樣本進(jìn)行模擬壓縮試驗(yàn)，計(jì)算樣本在z軸方向上的彈性模量，作為樣本的真值。在有限元仿真時(shí)，將骨骼材料視為線(xiàn)彈性、各向同性的均質(zhì)材料[12]，定義其彈性模量為18 GPa，泊松比為0.3[13]。骨小梁有限元仿真受力云圖，如圖4所示。

圖4 骨小梁樣本的受力云圖Fig.4 Stress results of trabecular bone samples

有限元分析得到樣本真值后，以24:1的比例將數(shù)據(jù)集隨機(jī)劃分為訓(xùn)練集（239 864）和測(cè)試集（10 000），隨機(jī)數(shù)種子設(shè)置為42。

3 結(jié)果

實(shí)驗(yàn)平臺(tái)配置如下：CPU為Intel(R) Xeon(R)E5-2630 v4 @ 2.20 GHz(×2)；顯卡為NVIDIA Tesla P100(×3)，顯存容量16 GB；操作系統(tǒng)為Windows Server 2012 R2 Standard，配置了CUDA 9.0和cuDNN 7.1.4；深度學(xué)習(xí)框架為Keras 2.2.4，以TensorFlow 1.12作為后端。

3.1 模型訓(xùn)練

使用HE[14]的初始化方法對(duì)模型的權(quán)重和偏差進(jìn)行初始化。為了防止過(guò)擬合，采用了L2正則化方法，正則化系數(shù)定義為10-4。模型以均方差作為損失函數(shù)，使用RADAM[15]方法進(jìn)行訓(xùn)練，初始學(xué)習(xí)率為10-3。受限于顯卡的數(shù)量和顯存容量，設(shè)置批量大小(Batch Size)為128×3。訓(xùn)練過(guò)程中，通過(guò)ModelCheckpoint、ReduceLROnPlateau、EarlyStopping等回調(diào)函數(shù)對(duì)測(cè)試誤差進(jìn)行實(shí)時(shí)監(jiān)測(cè)，若測(cè)試誤差比之前有所改進(jìn)，則將模型參數(shù)以.h5文件的形式保存到指定路徑；若測(cè)試誤差連續(xù)5輪沒(méi)有降低，則將學(xué)習(xí)率減少為原來(lái)的一半，以改善訓(xùn)練情況；若測(cè)試誤差連續(xù)15輪沒(méi)有降低，則終止訓(xùn)練，防止模型過(guò)擬合。經(jīng)過(guò)124輪的訓(xùn)練，模型的測(cè)試誤差停止下降。模型的訓(xùn)練誤差和測(cè)試誤差變化曲線(xiàn)，如圖5所示。觀察可知，訓(xùn)練誤差和測(cè)試誤差接近，模型沒(méi)有出現(xiàn)過(guò)擬合。

圖5 訓(xùn)練誤差和測(cè)試誤差變化曲線(xiàn)圖Fig.5 Curves of training error and test error

3.2 模型評(píng)估

SE-DenseVoxNet訓(xùn)練結(jié)束后，通過(guò)測(cè)試集對(duì)模型進(jìn)行交叉驗(yàn)證。針對(duì)骨小梁模量預(yù)測(cè)這一回歸問(wèn)題，使用均方根誤差(RMSE)、平均絕對(duì)誤差(MAE)、平均相對(duì)誤差(MAPE)和相關(guān)系數(shù)這4項(xiàng)指標(biāo)評(píng)估模型的表現(xiàn)。為了更好地對(duì)比模型的預(yù)測(cè)結(jié)果，實(shí)現(xiàn)了CIUSDEL等[8]的網(wǎng)絡(luò)模型，并使用上述4項(xiàng)指標(biāo)對(duì)其預(yù)測(cè)結(jié)果的準(zhǔn)確性進(jìn)行了評(píng)估，結(jié)果如表2所示。

表2 不同模型的評(píng)估結(jié)果對(duì)照表Tab.2 Comparison of evaluation results of different models

由表2可知，SE-DenseVoxNet的預(yù)測(cè)結(jié)果更準(zhǔn)確，平均相對(duì)誤差僅為1.779%，相比于CIUSDEL等的方法有明顯的提升。為了進(jìn)一步評(píng)估SE-DenseVoxNet的性能，對(duì)SEDenseVoxNet的預(yù)測(cè)結(jié)果及測(cè)試誤差進(jìn)行統(tǒng)計(jì)分析。SE-DenseVoxNet對(duì)測(cè)試集樣本的預(yù)測(cè)結(jié)果如圖6所示。其中橫坐標(biāo)為樣本真值，縱坐標(biāo)為網(wǎng)絡(luò)預(yù)測(cè)值，直線(xiàn)方程為y=x。觀察可知，樣本點(diǎn)分布在直線(xiàn)y=x附近，表明樣本真值和網(wǎng)絡(luò)預(yù)測(cè)值基本滿(mǎn)足線(xiàn)性關(guān)系。

圖6 真值-預(yù)測(cè)值散點(diǎn)圖Fig.6 Scatter plot of true values and predictions

使用Bland-Altman[16]方法對(duì)樣本真值和預(yù)測(cè)結(jié)果之間的一致性進(jìn)行分析。經(jīng)統(tǒng)計(jì)，樣本真值和預(yù)測(cè)結(jié)果的差值的平均數(shù)為1.468×10-4，標(biāo)準(zhǔn)差Sd為0.024，95%一致性界限為-4.703×10-2～4.732×10-2。根據(jù)統(tǒng)計(jì)結(jié)果，繪制了Bland-Altman散點(diǎn)圖。如圖7所示，散點(diǎn)圖的橫坐標(biāo)為樣本真值與預(yù)測(cè)結(jié)果的平均數(shù)，縱坐標(biāo)為樣本真值與預(yù)測(cè)結(jié)果的差值，中間虛線(xiàn)表示差值的平均數(shù)，上下兩條虛線(xiàn)表示差值的95%一致性界限。從圖7中可以看出，5.63%的樣本點(diǎn)在95%一致性界限以外；在95%一致性界限內(nèi)，差值的絕對(duì)值最大為0.047（白色標(biāo)記點(diǎn)），差值的平均數(shù)接近于0，說(shuō)明樣本真值和預(yù)測(cè)結(jié)果之間的偏倚很小，具有較好的一致性。

圖7 Bland-Altman散點(diǎn)圖Fig.7 Bland-Altman scatter plot

4 討論

提出了一種基于SE-DenseVoxNet的骨小梁模量預(yù)測(cè)方法，并通過(guò)骨小梁三維圖像數(shù)據(jù)集對(duì)網(wǎng)絡(luò)進(jìn)行了交叉驗(yàn)證。實(shí)驗(yàn)證明，本方法的預(yù)測(cè)結(jié)果與樣本真值之間的誤差和偏倚很小，具有較好的一致性。相比于CIUSDEL等的方法，基于SE-DenseVoxNet的模量預(yù)測(cè)方法可以更準(zhǔn)確地預(yù)測(cè)樣本的生物力學(xué)性質(zhì)，輔助醫(yī)生診斷病人的骨骼健康狀況。在接下來(lái)的研究中，會(huì)嘗試獲取更多的骨骼樣本，擴(kuò)充數(shù)據(jù)集，搭建更復(fù)雜的網(wǎng)絡(luò)模型，進(jìn)一步提升模型的預(yù)測(cè)精度，減少模型誤差。