Segmentation of retinal OCT images using a random forest classifier

Segmentation of retinal OCT images using a random forest classifier

使用隨機森林分類器對視網膜OCT影像進行分割

抽象

光學相干斷層掃描（OCT）已成為診斷視網膜異常的最常用工具之一。視網膜形態和層厚均可提供重要資訊，以幫助對這些異常進行鑑別診斷。自動分段方法對於提供這些厚度測量至關重要，因為考慮到每次OCT掃描中的資料量很大，每層的手動描繪都很麻煩。在這項工作中，我們提出了一種使用隨機森林分類器進行視網膜層分割的新方法。從OCT資料中總共提取了七個特徵，並用於同時對九個層邊界進行分類。利用隨機森林的概率性質，提取每個邊界的概率圖，並將其用於幫助完善分類。我們能夠準確地分割八個視網膜層，其平均Dice係數為0.79±0.13，並且層邊界的平均絕對誤差為1.21±1.45畫素。
關鍵詞
OCT; 視網膜層分割 隨機森林分類

1.引言

光學相干斷層掃描（OCT）是一種新興技術，在眼科學1中還有許多其他用途。 OCT的一種新興應用是評估多發性硬化症（MS）患者。 已顯示患有MS的受試者表現出各種視網膜層的顯著變薄。2自動分割方法對於估計這些視網膜厚度變得尤為重要，因為每層的手動分割是一項耗時的任務。 迄今為止，在文獻中已經提出了各種各樣的自動視網膜分割方法。 許多作者提出了自動視網膜層分割方法，利用了多種技術，包括活動輪廓3，基於梯度的資訊4，5，統計模型6，基於圖的方法7-9以及畫素分類10。
基於圖形的OCT分割方法越來越受歡迎。 視網膜的分層結構非常適合這種問題的形成。 這些方法還被證明是高效，簡單且具有良好的最優性。 Chiu等人[8]從每張2D OCT影像中構造出一個圖形，並使用最短路徑演算法分割多層。 Garvin等人7構建了4D圖，以對3D OCT體積和所有圖層進行建模。 使用新穎的圖形切割公式解決了該問題。 這些基於圖的方法的問題之一是在圖上找到適當的成本函式來區分每一層。 Antony等[9]在先前的工作7的基礎上，通過使用近鄰分類器將影像特徵分組到不同的區域中，然後將其用於圖形權重中。
使用畫素分類進行分割的另一種方法來自Vermeer等人[10]。這項工作使用一組18個梯度和強度特徵，將其投影到更高維度的空間中，為每個視網膜層建立支援向量機分類器。 然後使用級別設定方法對獨立分類的畫素進行正則化以建立平滑表面。 對每個層使用單個分類器，對每個表面使用單獨的正則化方法，會導致效率特別低下。 另外，僅給出了五層介面的分割結果。
在這項工作中，我們提出使用隨機森林分類器 11自動了解黃斑OCT影像中適當的邊界在哪裡。 總共檢測到九層邊界，導致八層分割。 這種型別的分類器在不同的影像分割任務12、13中表現出良好的效能，並且由於易於使用，需要調整的引數數量少以及能夠同時對多個標籤進行分類而廣受歡迎。 隨機森林的一個特別好的功能是，它以概率的形式為每個標籤提供了軟分類。 在對影像進行分類時，可以生成概率圖，然後將其細化為硬分割。 在視網膜分割的情況下，我們在邊界而不是層上進行分類，然後可以利用層的有序結構，通過為每個邊界每次A掃描僅選擇一個畫素來將概率細化為分割。

2.方法

提出的工作目的是在黃斑OCT掃描中分割八個視網膜層（圖1）。從視網膜的頂部到底部列出要分割的層是：視網膜神經纖維層（RNFL），神經節細胞層和內部叢狀層（GCL + IPL），內部核層（INL），外部叢狀層（OPL），核外層（ONL），內部部分（IS），外部部分（OS）和視網膜色素上皮（RPE）。這八個層有九個關聯邊界，請注意，內部限制膜（ILM）用於表示玻璃體RNFL邊界，外部限制膜（ELM）用於ONL-IS邊界和Bruch膜（BM）用於RPE脈絡膜邊界。單個2D OCT影像稱為B掃描影像，其中A掃描是B掃描中的單個畫素列。多個B掃描在平面方向上組合以形成3D體積。我們還將x和y方向分別表示為B掃描中的從左到右和自上而下的方向，將z方向表示為穿過平面的方向。
在分割視網膜層之前，會自動生成視網膜蒙版，以粗略定義我們希望在其中找到這些層的感興趣區域。 視網膜掩模的計算還使我們能夠將每個影像展平到BM。 為了分割視網膜層，從OCT資料中提取了七個特徵，並用於訓練隨機森林分類器以找到每個邊界。 分類器產生每個畫素屬於每個邊界的概率。 然後對這些概率進行細化，以估計每層的最終分割。

2.1視網膜檢測和展平

在分割視網膜層之前，我們會生成一個粗糙的視網膜遮罩，以指示哪些畫素在視網膜的內部和外部。 圖2（b）示出了圖2（a）中的OCT影像的示例性視網膜掩模，其中白色和黑色分別表示視網膜內部和外部的區域。 視網膜掩模的計算需要ILM和BM邊界的估計。 由於這是一個預處理步驟，因此需要快速計算。 此外，由於這些邊界將在以後的階段中進行完善，因此只需大致定位即可。
為了計算視網膜遮罩，首先使用高斯濾波器（σx=σy= 10）平滑體積中的每個B掃描影像。沿每次A掃描檢視，假定具有最大正梯度值的畫素為ILM或IS-OS邊界。在距先前找到的最大值至少25個畫素處具有最小正梯度值的最大正梯度值的畫素被視為第二邊界。給定這兩個畫素，最接近影像頂部的一個作為ILM。然後，沿著每次A掃描，將BM作為IS-OS以下的最大負梯度值。由於這些估計的ILM和BM表面可能包含虛假的跳躍和不連續性（例如，由於血管偽影所致），因此我們刪除並用最近的點填充外圍邊界點。外圍點是那些距離其各自的10×10中值濾波表面超過15個畫素的點。最後，使用高斯核對兩個表面進行平滑處理（對於ILM，σx，z = {10，0.75}，對於BM，σx，z = {20，2}）。然後，視網膜蒙版體積包含估計的ILM和BM表面之間的所有畫素。圖2（c）顯示了一個B掃描影像，其中非視網膜區域被遮蓋，表明該視網膜遮罩僅粗略地定位了頂部和底部邊界。
給定視網膜掩模，然後通過平移每個單獨的A掃描將OCT資料平整到BM，以使該邊界完全平坦。 雙線性插值用於平移。 平坦化過程去除了每個影像中的大部分曲率，將所有視網膜影像放置在跨物件的公共空間中，這是文獻中常見的一個步驟。5，7圖10和圖10是OCT影像及其產生的平坦化影像的示例。 2（a）和2（d）。

2.2隨機森林分類器

作為分割的第一步，訓練一個隨機森林分類器11以找到每個圖層的邊界畫素。 僅使用一個分類器來學習所有邊界，每個畫素被分配一個0到9的標籤，具體取決於該畫素是否位於九個邊界之一（標籤1–9）或不位於（九個標籤）0上。 森林中使用了五十棵樹，並且在每個節點上隨機選擇了兩個特徵。 為了訓練分類器，使用了七個特徵（如圖3所示）。 前兩個特徵是位置特徵，後五個特徵是區域性影像特徵。 這些特徵是在平坦的OCT影像中的每個畫素處計算的。
第一個特徵是利用Sec21中計算的視網膜蒙版為視網膜所特有的。 它是從視網膜蒙版的底部到頂部的每個畫素的分數距離，取值介於0和1之間。視網膜蒙版上方或下方的畫素被分配為-1的值，以幫助分類器更好地區分視網膜與完全背景。 下一個使用的特徵是影像中每個畫素的x位置（在我們的資料中採用1到1024的整數值）。 此功能的目的是幫助瞭解整個影像上層厚度的變化。 它與視網膜掩模距離互補工作，以區分高曲率中央凹區域，由於掃描協議，中央凹區域始終朝著影像中心出現。 這兩個位置特徵在圖3（b）和3（c）中示出。
接下來，通過獲取每個畫素上方和下方3×3鄰域的平均值之差，計算出四個多尺度，類似梯度的特徵。 鄰域位於上方和下方的2、4、8和16畫素處，這有助於解釋周圍圖層的強度差異。 這些特徵的示例如圖3（d-g）所示。 此功能旨在瞭解附近有哪些圖層（即上方的暗層和下方的亮層）。 最終特徵是每個畫素上結構張量的最大特徵值的大小。14張量積用sigma值0.5平滑。 每個畫素的第一特徵值在各向異性區域較大，因此特別擅長在視網膜底部找到較薄的層（圖3（h））。
接下來，通過獲取每個畫素上方和下方3×3鄰域的平均值之差，計算出四個多尺度，類似梯度的特徵。 鄰域位於上方和下方的2、4、8和16畫素處，這有助於解釋周圍圖層的強度差異。 這些特徵的示例如圖3（d-g）所示。 此功能旨在瞭解附近有哪些圖層（即上方的暗層和下方的亮層）。 最終特徵是每個畫素上結構張量的最大特徵值的大小。14張量積用sigma值0.5平滑。 每個畫素的第一特徵值在各向異性區域較大，因此特別擅長在視網膜底部找到較薄的層（圖3（h））。

2.3最終細分

隨機森林分類器為每個畫素輸出每個標籤的概率，該概率是根據投票給該標籤的隨機森林中決策樹的數量計算得出的（見圖4）。 由於我們希望每個A掃描僅給每個邊界標籤分配一次，因此在每個畫素處以最大概率獲取標籤效果不佳。 相反，對於每個類別的概率影像，我們每次A掃描都仔細選擇一個畫素位於其關聯的邊界上。
為了從圖層的概率影像中提取單個邊界，首先將高斯平滑（σ= 15）應用於影像。沿每次A掃描具有最大概率的畫素被選為該類別的初始邊界畫素。如圖4所示，由於這些邊界通常很嘈雜，因此將概率低於閾值0.05的畫素刪除，然後使用具有較大平滑引數的三次平滑樣條對其餘點進行平滑處理。可以對每一層重複此過程以估計所有邊界。最後，在適當的邊界之間填充圖層標籤。為了確保各層的順序正確，從頂層和底層開始，一次填充一層。在新增每個圖層時，將施加硬約束，以使新圖層無法相交併穿過先前擬合的圖層。如果曲面位於先前擬合曲面的外部，我們將調整這些點以在相鄰邊界上方或下方放置一個畫素。仔細選擇如何新增圖層很重要，因此首先新增具有較高置信度的圖層。

3.結果

在這項研究中使用了十個手動繪製的OCT卷。 所有手動描述均由單個評分者使用內部開發的軟體完成。 提供知情同意後，患者在Spectralis OCT掃描器（德國海德堡海德堡工程公司）上進行光譜域視網膜成像。 獲得由49次B掃描組成的黃斑光柵掃描（20°×20°）。 該掃描器利用自動實時（ART）技術，通過採集和平均同一位置的多幅影像並進行平均以降低噪聲，從而提高影像質量。 本研究使用ART≥12（平均掃描數）且訊號質量至少為20 dB的掃描。 影像體積每個A掃描有496個畫素，每個B掃描有1024個A掃描。 還要注意，原始頻譜強度資料的第四次根是在所有處理之前應用的。
五次OCT卷中的一小部分資料用於訓練隨機森林分類器。在每個訓練量中，選擇了三個B掃描作為代表掃描。三個掃描中的一個掃描包括中央凹中央區域，而其他兩個掃描則從中央凹的每一側進行。最後，除了在每個邊界之間找到相等數量的隨機背景點之外，每個邊界上的所有點（每次B掃描每層1024個點）都用於訓練。其餘5個OCT卷用於探討演算法的執行情況。使用兩種方法來分析演算法的效能。首先，我們看一下每個估計的層邊界與地面真實手動分割之間的平均絕對誤差（MAE）。使用的第二個量度是Dice係數，它是估計的層分割與真實的層分割之間的重疊的量度。對於每個層，將其計算為，其中Xk和Yk分別是在地面真實層和估計層中標記為k層的點集。

表1和表2顯示了我們演算法在所有物件和B掃描中的平均MAE和Dice係數結果。根據表1，ILM的邊界誤差最低。可以預期該結果，因為該邊界是背景與視網膜相交併具有較大漸變的位置。另外，參見圖2（c），視網膜掩模非常緊密地定位ILM，從而使視網膜掩模功能在該邊界上非常好。再次由於較大的可識別梯度，ELM，IS-OS和BM邊界也顯示出較低的誤差。 OS-RPE邊界被證明是最難估計的，因為該邊界傾向於消失並重新出現在某些影像中。從Dice係數來看，所有層均顯示出相當不錯的結果，而IS，OS和RPE層則存在更多困難。 Dice係數最大的圖層RNFL，GCL + IPL和ONL比許多其他圖層都寬，並且在邊界處也具有較大的漸變。由於IS，OS和RPE的層較薄，因此它們的效能最差。這些層中的一個畫素邊界誤差轉化為比在較大層中具有類似邊界誤差的更低的Dice係數。
圖5顯示了我們的分割方法在三個B掃描上的結果。 第一行顯示輸入的B掃描，第二行顯示手動分割，最後一行顯示自動分割的結果。 一個明顯的困難領域是在包含血管偽影的區域。 如圖4所示，在概率影像中，帶有血管的區域要麼顯示為假陽性，要麼接收到較低的概率值。 誤報會導致邊界不正確，然後在最後的細分步驟中將邊界平滑。 但是，這種平滑效果不足以消除錯誤。 當邊界的可能性較低時，通常會有足夠的響應來使邊界合理匹配。 如果響應太小（低於0.05閾值），則由於表面平滑變化，因此平滑處理可以非常準確地對缺失區域進行插值。
最後，我們快速檢視演算法的效能。 所有實驗均在執行Windows 7的計算機上執行，該計算機使用MATLAB R2011b在1.73 GHz頻率下執行四核處理器。 隨機森林實現可線上線上獲取。15此實現利用預編譯的mex檔案來提高效率。 所有其他處理都使用本機MATLAB函式。 該演算法不同方面的效能如表3所示。訓練分類器耗時約90秒。 如圖所示，分割一個卷的總時間約為117秒，包括所有步驟。 此時間對應於每次B掃描2.4秒。 使用優化的程式碼可能會進一步提高效能。 例如，已在GPU上實現了隨機森林的有效實現。16
4。結論
我們展示了僅使用七個特徵就可以使用隨機森林對視網膜層進行分類的能力。 為了為每一層生成最終的分割，需要對隨機森林分類器輸出的邊界概率影像進行後處理。 該演算法被證明是快速而準確的，在僅對五個體積的一小部分進行訓練後，就可以成功分割多個OCT體積。 這項工作有許多方面可以改進。 該演算法本質上是2D的，在單獨的B掃描上執行，而不是使用來自完整3D卷的資訊。 使用全部體積可能有助於解決血管偽影的問題。 還必須執行進一步的驗證，比較評估者之間和評估者內部的細分差異。

圖1手動裁剪的OCT影像 九個分段的邊界從上到下分別是：ILM，RNFL-GCL，IPL-INL，INL-OPL，OPL-ONL，ELM，IS-OS，OS-RPE和BM。

圖2.（a，b）一張OCT影像和計算出的視網膜掩模。 （c）遮蓋了非視網膜畫素的OCT影像，顯示出粗略的準確性。 （d）OCT影像展平至底部邊界。

圖3. OCT影像（a）及其相關特徵：（b）視網膜掩模距離，（c）x位置，（d–g）4個遞增比例的梯度特徵，以及（h）結構張量的最大特徵值。

圖4.從（a）的影像上的隨機森林分類器生成的所有九層的OCT影像（a）和概率圖（bj）。 概率介於0（黑色）和1（白色）之間。

圖5.三個B掃描的分割結果，每列表示一個不同的B掃描。 頂部，中間和底部各行分別是輸入的OCT影像，手動分割和自動分割。

表1 每個邊界的平均絕對誤差（標準偏差）值。 值以畫素為單位。

表2估計的層分割和真實層分割之間的Dice係數值。

表3演算法各個方面的執行時間