基于旋轉不變非線性局部模糊編碼的膝骨關節炎輔助診斷

任會峰，伏建雄，鄢鋒，李富

（1.無錫學院物聯網工程學院，江蘇無錫 214105；2.湖南中自工業智能研究院有限公司，湖南長沙 410000；3.長沙有色冶金設計研究院有限公司，湖南長沙 410011）

0 引言

膝骨關節炎（Knee Osteoarthritis，KOA）是世界上致殘率第四高的疾病，其病因和發病機制尚不明確，性別、職業、年齡、遺傳、運動過度、體重超重等多種因素都可以誘發或加重KOA的進展。控制膝骨關節炎的最佳策略是早發現、早診斷、早治療。研究表明，磁共振成像（Magnetic Resonance Imaging，MRI）對眾多疾病的診斷表現出良好的敏感性、特異性和可重復性。毛家驥等依據MRI加權像分析手腕部關節炎表現，確定有關病證分布及發生規律；王紫儀等探討磁共振成像加權像高信號比值的病因推測價值及其與聽力預后的關系，證明信號強度越高，提示出血可能性越大，聽力預后越差。尤其在器官磁共振評分系統中，通過定量測量T2值，可以量化結構成分和生物力學指標變化，相較于關節鏡等檢查方法，磁共振檢查方式準確率更高、漏診率更低，組間差異顯著。因此，T2加權像表現出與膝骨關節炎較強的相關性。依據磁共振加權像推測膝骨關節炎病況的傳統方式主要依賴于醫生的主觀經驗，對醫生接診的病例類型數目、積累的醫療診斷經驗等要求較高，受專業水平、情緒波動等多種因素影響難免存在誤判，容易導致醫療糾紛，而且診斷效率降低。智慧醫療技術的發展為眾多疾病提供了客觀的診斷依據和高效的治療手段。相關領域廣泛采用圖像處理技術輔助疾病診斷，受疾病醫學圖像樣本數量少，尤其是難以獲得數據完備的樣本，以及圖像標注缺乏領域知識等因素影響，一些基于大量數據樣本的圖像處理技術受到限制。

紋理是一種無處不在的視覺特征，已廣泛應用于多種疾病的輔助診斷。劉欣等應用Logistic回歸分析行多特征聯合回歸分析相關紋理特征ROC曲線的AUC，論證了MR T2WI圖像紋理對于鑒別類風濕關節炎和痛風性關節炎具有一定價值。研究表明，紋理分析能充分挖掘醫學影像中肉眼無法觀察到的隱藏信息，能反映組織的潛在病理及生理學特征，在鑒別診斷方面潛力巨大。基于圖像紋理的計算機輔助診斷技術用于膝骨關節炎檢測，可以減少對診斷的主觀性影響，降低對操作者熟練程度的依賴性，避免因疲勞等因素造成誤診。紋理特征提取為紋理分析提供數據基礎，從實用角度講，要求紋理特征數量少、區分度高、魯棒性好且易于理解。目前主流的紋理提取方法有基于統計的方法、基于建模的方法和基于頻域變換的方法。然而，由于運動模糊、斑點噪聲等偽影造成膝關節T2加權像的不確定性或不完備性，目前流行的灰度共生矩陣、分形計算、Gabor濾波等方法性能不佳。而紋理直方圖計算簡單、魯棒性強、易于理解，是描述局部區域紋理的有效方法。湯衛霞等利用MR圖像的紋理直方圖特征揭示腫瘤異質性，表明MR圖像紋理直方圖在腫瘤的早期診斷、精準治療和預后評估等方面具有巨大潛力。但是，紋理直方圖的傳統局部硬編碼基于鄰域像素和種子像素的二值比較提取紋理特征，其前提是不同紋理單元的紋理單元向量不同，這一前提假設與實際不符，而且局部硬編碼忽略了像素分布的不確定性和灰度差導致紋理單元數量減少隨之圖像區分能力下降。在此基礎上，Barcelo等提出采用局部模糊編碼研究骨髓顯微圖像中白細胞和紅細胞圖像變化和內在空間結構敏感性，其采用線性隸屬度函數，提高了圖像分類準確率，但是利用線性隸屬度函數擬合像素灰度差的包絡線存在較大偏差，而且得到6 561維的高維特征向量描述局部模糊編碼，特征維數偏高、冗余噪聲較大且對旋轉圖像分類效果不佳。

本文在分析傳統局部硬編碼的基礎上，采用指數函數代替線性隸屬度函數擬合鄰域像素與中間像素灰度差的包絡線，通過局部模糊編碼頻率直方圖統計描述局部單元的灰度變化，借鑒串移位操作將旋轉后局部模糊編碼相同的局部單元歸為一類，構建一種旋轉不變的非線性局部模糊編碼紋理描述方法，從眾多局部單元中選取最具代表性的局部單位向量描述圖像紋理，從而大幅降低紋理特征的維數。在此基礎上，通過構建一個局部模糊編碼描述符進一步壓縮紋理特征的維數，并用于紋理圖像分類和膝骨關節炎輔助診斷。

1 傳統局部硬編碼

傳統局部硬編碼（Local Grid Coding，LGC）是基于圖像塊鄰域像素灰度值大小關系的紋理描述算法，其基本思想是：任意圖像

I

均由多個3×3共9個像素大小的局部單元組成，將局部單元中間像素

V

灰度值記為

V

，周圍的8個領域像素

V

到

V

灰度值記為

V

，

i

=1…8。假設

V

為局部單元左上角像素，從

V

開始按照順時針方向將每個鄰域像素與中間像素的灰度值進行比較，得到一個8維的局部單元向量（Local Unit Vector，LUV），

LUV

=（

E

，…，

E

），其中：

根據式（1）可知，局部單元可能有3～3種

LUV

，不同局部單元的

LUV

數用局部單元值

N

區分，其中：

所有局部單元值的頻率直方圖即為傳統局部硬編碼?？梢姡植坑簿幋a需要高達6 561維的單元向量描述，而且，從局部硬編碼的計算過程可知，該算法忽略了像素分布的不確定性和像素灰度值差異程度。為了描述鄰域像素灰度值差異程度，文獻［24］提出線性局部模糊編碼（Linear Local Fuzzy Coding，LLFC），但該算法采用的模糊隸屬度函數為線性，利用線性隸屬度函數擬合鄰域像素灰度差異大小的包絡線存在較大偏差，因此本文提出非線性局部模糊編碼。

2 非線性局部模糊編碼

為便于描述，將局部硬編碼中

E

隸屬于0、1、2的第

i

個鄰域像素的隸屬度值記為

μ

（

d

）、

μ

（

d

）和

μ

（

d

），構成局部模糊單元向量

FUV

。

由式（3）可知，鄰域像素

V

和中間像素

V

之間的灰度差

d

越小，鄰域像素和中間像素灰度值越相似，反之亦然，則圖像

I

鄰域像素與中間像素灰度差直方圖的包絡線如圖1中實線所示，采用指數函數對該包絡線進行擬合，圖1中虛線為指數擬合曲線。

Fig.1 Envelope and fitting curve圖1 包絡線與擬合曲線

從圖1可知，|

d

|＞40的概率小于0.1，故可設置非線性模糊隸屬度函數如下：

由式（1）-式（6）所示，對任意圖像

I

，從圖像

I

的第2行第2列像素作為中間像素開始到倒數第2行第2列像素為止，以3*3的滑動窗口掃描整幅圖像以計算所有局部單元值出現的頻率，圖像

I

的非線性局部模糊編碼就是所有局部單元值的頻率直方圖。可見，非線性局部模糊編碼很好地消除了局部硬編碼采用絕對差值的缺陷，但仍然需要很高的維數向量以描述局部模糊編碼，最高可達6 561維，而且該算法依然沒有考慮圖像角度旋轉問題。文獻［25］指出，對稱不變的紋理模式應該視為同一紋理模式，對一幅圖像旋轉一定角度后，圖像的紋理并沒有發生變化，但是按照非線性局部模糊編碼的計算過程卻會出現完全不一樣的非線性局部模糊編碼，因此傳統局部硬編碼以及非線性局部模糊編碼均不具備旋轉不變性。此外，高達6 561維的向量不僅導致較大的計算量，而且會為圖像的解讀增加干擾噪聲，影響機器的判斷。為此，對非線性局部模糊編碼加以改進，得到旋轉不變非線性局部模糊編碼（Rotation Invariant Nonlinear Local Fuzzy Coding，RLFC）。

3 旋轉不變非線性局部模糊編碼

圖2為任意圖像中都可能出現的12種局部硬編碼的局部單元，大致分為圖2（a）和圖2（b）兩大類。借鑒字符串中的移位操作可知，在順時針右移位之后，圖2（a）中的局部單元向量可以轉換為下一個局部單元向量。在圖2（b）中也有相同的結構。

Fig.2 Examples of local unit圖2 局部單元示例

在此，將右移后局部單元向量相同的局部單元歸為一類，則局部單元維數從6 561維縮減到834維，大幅削減描述圖像紋理的特征維數，而且具備旋轉不變性，描述圖像紋理的能力增強，但描述圖像紋理的向量維數依然較高。

因此，

f/f

比值越小，圖像越粗糙。當

f

＜0.01*

m*n

時，圖像極度粗糙。

f+f

值越大，圖像越粗糙；

f/f

比值越小，圖像越平滑，當

f

＜0.01m*n，圖像極度平滑。

f+f

值越大，圖像越平滑，其中，

m*n

為圖像的分辨率。記

f

=

f

+

f

，

f

=

f

+

f

，可知，對于任意圖像，若

f

越大，

f

越小，則圖像越粗糙；

f

越大，

f

越小，圖像越平滑。從而，可用一個2維向量

X=

（

f

，

f

）表示局部模糊編碼直方圖的復雜度，大幅降低了描述圖像紋理的特征維數，且削減了干擾信息。

4 實驗結果及分析

4.1 標準圖像集測試

為了評價本文所提方法的有效性，將所提算法RLFC與傳統局部硬編碼LGC、線性局部模糊編碼LLFC、非線性局部模糊編碼NLFC進行比較，從Kylberg紋理圖像庫和Brodatz紋理圖像庫中抽取測試圖像集，這兩個圖像庫常用于評價和比較紋理度量。其中，Kylberg圖像庫包括28個類，每個類有160個圖像，每幅圖*像的大小為576×576像素，圖像數據集可從http：//www.cb.uu.se/～gustaf/texture下載，從每個類中隨機抽取一幅圖像均分為4個大小為288×288像素的子圖像，將每一幅子圖像以10°為間隔，旋轉10°到360°，從旋轉后的子圖像中心提取大小為128×128像素的孫圖像，共得到28×4×36=4 036幅圖像，作為本次測試的Kylberg圖像集；Brodatz圖像庫可從http：//www.ee.oulu.fi/mvg/page下載，包括112類紋理圖像，從每一類中隨機抽取一張紋理圖像以10°為間隔，旋轉10°到360°，從旋轉后的紋理圖像中心提取大小為128×128像素的子圖像，共得到112×36=4 032幅圖像，作為本次測試的Brodatz圖像集。采用式（8）基于

L

-距離的

K

-means聚類：

其中，

X

、

X

是兩幅圖像的局部模糊編碼描述向量，

x

、

x

是

X

和

X

的分量。將聚類結果與原始標簽一一對比，正確歸類的樣本數與總樣本數之比定義為精度，以此作為評價指標。隨機初始化聚類中心，并分別執行10次聚類計算。表1為在構建的Kylberg圖像集和Brodatz圖像集上分別采用傳統局部硬編碼LGC、線性局部模糊編碼LLFC、非線性局部模糊編碼NLFC和本文提出的旋轉不變非線性局部模糊編碼RLFC聚類精度對比。

從表1可以看出，局部硬編碼可以在一定程度上反映鄰域像素和種子像素的關系，但灰度值的硬劃分導致聚類精度較低。此外，對于不同的初始聚類中心，該算法表現出較弱的魯棒性。然而，無論是否對紋理單元向量進行三進制編碼，局部模糊編碼、旋轉不變非線性局部模糊編碼都能較大程度地提高準確率，在經過旋轉預處理所構建的圖像數據集上，本文提出的方法相比局部模糊編碼表現出更高的準確率。10次測試中，最大準確率分別提高9%和1%，最小準確率分別提高6%和3%，平均準確率分別提高14%和7%。這是因為本文提出的方法雖然減少了紋理單元數，但通過局部單元平移，提高了圖像的旋轉不變性和對不同初始聚類中心的魯棒性。

Table 1 Comparison of Clustering accuracy of different algorithm on texture ionage dataset表1 不同算法在紋理圖像集上的聚類精度對比

為了進一步體現所提方法的優越性，在上一步構建的Kylberg和Brodatz圖像集上添加均值為0、信噪比（Signal to Noise Ratio，SNR）由20變化到10的高斯白噪聲，以評價各算法抗噪聲干擾能力，依然采用式（8）的K-means聚類，以準確率作為評價指標，圖3（a）和圖3（b）分別為各算法在增加白噪聲后的Kylberg和Brodatz圖像集的實驗結果。

從圖3可知，隨著SNR的減小，噪聲水平逐漸增加，各算法在兩個圖像集上的變化趨勢相同，聚類的準確率都呈下降趨勢。其中，傳統局部硬編碼LGC的聚類精度下降最快，表明該方法的抗噪聲能力最差；線性局部模糊編碼和非線性局部模糊編碼在兩個數據集上隨著噪聲水平的增加，聚類精度略有下降，但總體變化較小，所提方法RLFC在不同高斯白噪聲影響下的性能均優于其它方法，而且在增加高斯白噪聲的兩個紋理圖像集上表現都比較穩定。綜合表1和圖3可見，所提方法無論在聚類精度還是在抗噪聲干擾能力上都較好。

Fig.3 Experimental results after adding Gaussian white noise圖3 增加高斯白噪聲后的實驗結果

4.2 膝骨關節炎輔助診斷

圖4依次為KOA逐漸嚴重的4幅T2加權像，與圖像對應的病癥都具有明確的病理診斷結果，而且相關患者在MRI檢查前未接受穿刺取樣及手術治療，將圖像進行灰度均衡化處理，以降低不同設備不同參數的圖像對比度和亮度對圖像灰度值的影響，并由具有多年MRI診斷經驗的醫師進行閱片分析及標注KOA分期，對應KOA的初期、中期、后期和晚期，4幅圖像中存在不同程度的關節軟骨損傷，但對病變嚴重程度肉眼難以區分。

Fig.4 Images of KOA圖4 KOA樣本圖像

圖5為基于本文提出方法的4幅KOA圖像的局部模糊編碼直方圖。從圖5可以看出，不同的KOA圖像在局部模糊編碼描述符上有比較明顯的差異。為了進一步直觀展示不同樣本圖像之間的紋理差異，按照式（7）對各圖像分別計算4個典型局部單元類及出現的頻率。

Fig.5 Local fuzzy coding histogram of T2 weighted image圖5 T2加權像的局部模糊編碼直方圖

表2是各幅圖像局部模糊編碼直方圖復雜度的二維向量（

f

，

f

）對照表。從表2可見，KOA1的

f

最大，

f

最小，表明KOA1對應的圖像紋理最平滑，這是因為初期血管破壞較少，軟骨纖維化較輕。從初期到晚期，

f

逐漸減小，

f

逐漸增大，KOA4對應的圖像紋理最粗糙，這是因為晚期滑膜破壞超過全面積的75%，軟骨纖維化嚴重，甚至形成骨質增生骨贅，血管破壞嚴重到血管閉塞。

Table2 RLFC descriptor of KOA images表2 KOA樣本圖像的RLFC描述符

圖6為

f

/

f

比值和KOA嚴重程度之間的關系曲線。由圖6可知，從初期到晚期，

f

/

f

比值逐漸減小，進一步說明隨著膝骨關節炎的不斷惡化，T2加權像的紋理越來越粗糙。從曲線的斜率變化趨勢看，KOA發展的4個階段中，從初期到中期變化的過程相對緩慢，這是由于KOA屬于多因素導致的相應部位病理產物對人體刺激的綜合結果，這與初期人體生理組織相對較好，發生病變的組織相對較少有關，一般由于膝關節后側半膜肌腱的薄弱部位發生病變有關；從中期到后期惡化最為迅速，這是因為中期伴隨外側半月板前外緣游離，容易在扭轉時發生橫形破裂等不可逆的繼發損傷，而這種損傷會逐漸加??；后期一般發生于前側臏腱勞損，從后期到晚期發展也比較迅速，在后期發生髕骨軟化后膝關節內側半月板由于失去保護而加劇磨損形成KOA晚期。

Fig.6 Relationship between f T/f H ratio and KOA severity圖6 f T/f H比值與KOA嚴重程度的關系曲線

由表2和圖6可知，

f

越小，則

f

越大。KOA圖像的平滑度隨

f

/

f

的減小而增加，對應關節炎更為嚴重，測試結果與實際情況相吻合，表明所提出的方法對膝骨關節炎的分期和診斷有較好的輔助作用。

5 結語

采用MRI的紋理分析有助于術前評估膝骨關節炎的病理分期，是一種KOA輔助診斷的有益探索。局部模糊編碼摒棄了傳統紋理譜中的剛性比較，能夠更準確地描述圖像的局部特征。通過引入移位操作和三進制編碼，降低了紋理向量的維數，減少了計算量，且具有更好的旋轉不變性和抗噪性，在標準圖像集和膝骨關節炎T2加權像的實驗結果表明，提出的旋轉不變非線性局部模糊編碼的紋理分析方法是一種無創性KOA輔助分期手段，可以為KOA早期發現提供快速的客觀依據，有利于早期指導患者生活避免病程快速發展。

本研究為單中心的回顧性研究，雖然對獲取的實驗圖像進行了預處理等，一定程度降低了圖像數據的差異性，但僅僅基于MRI圖像紋理分析，并未結合臨床特征，如性別、年齡等，有待后續對影像數據進一步完善后結合更深入的影像組學方法與大數據關聯挖掘進行綜合分析。作為一個相對新興的領域，基于紋理分析的KOA輔助診斷處于起步階段，仍面臨很多挑戰，需加強實驗和數據處理技術以增加實驗樣本數，研究小樣本機器學習的輔助診斷技術和采用半自動或全自動方法實現病變位置的自動定位和病理解釋將是下一步工作的重點。