機器學習:U-NET ConvNet用於CT掃描分割

U-NET是基於卷積神經網絡(CNN),能夠產生視覺信息。通過計算訓練和明確定義的優化公式,可以獲得合理的結果(Dice Score約為0.9),同時在CT掃描上識別骨骼和肝臟,而無需對超參數進行進一步手動調整。進一步調整可能會帶來更好的結果。

介紹

計算機斷層掃描掃描產生的斷層掃描圖像允許我們從內部看到感興趣的物體而沒有開口。該技術廣泛用於幾種不同情況下的醫學診斷。

來自人體的斷層圖像通常顯示幾個器官,並且需要專家來分割感興趣的區域並解釋結果。

圖像分割是圖像處理的基本活動之一,通常是進一步分析的先前步驟。在CT掃描的情況下,醫生通常有興趣找到可以解釋患者症狀或狀況的異常,因此,他們必須能夠識別感興趣的部分並指出有關所研究結構的相關事實。

在這篇簡短的文章中,我們將通過卷積神經網絡探索醫學圖像分割。特別是,我們將研究U-NET結構並將其應用於真實的CT掃描數據並檢查其執行情況。

背景理論

人工神經網絡可以追溯到20世紀40年代,其中McCulloch和Pitts研究導致了一種神經元模型,該模型是從Rosenblatt的感知思想發展而來的。

1989年晚些時候,數學家Cybenko證明了具有sigmoid激活函數的人工神經網絡的通用逼近。換句話說,這基本上意味著具有有限數量的神經元的人工神經網絡能夠在緊湊的n維實數子集上產生連續函數的近似。

同樣在1989年,LeCun發表了他關於卷積神經網絡的著作。簡而言之,CNN使用卷積運算來計算神經網絡的權重,從而減少跨層的維數,同時保持重要的空間關係。

GPU技術的進步使基於CNN的技術在ImageNet等圖像識別競賽中達到了新的水平。

簡單來說,卷積是一個數學運算,它接受兩個函數f和g,併產生第三個h函數,描述一個函數的形狀如何影響另一個函數。

離散卷積是通過一個公式實現的。

機器學習:U-NET ConvNet用於CT掃描分割

U-Nets出現在2015年Ronneberger等人的文章中(https://arxiv.org/abs/1505.04597)。並且在2016年,Christ等人(https://arxiv.org/pdf/1610.02177.pdf)在CT掃描圖像上進行自動肝臟分割。關於U-Net的一個很好的想法是,它可以接收圖像作為輸入,生成另一個圖像作為輸出,這對於生成分割圖像非常方便。

機器學習:U-NET ConvNet用於CT掃描分割

2015年Ronneberger等人展示的U-Net

U-Net的思想是,通過訓練,網絡將能夠通過最小化與所需操作相關的成本函數來綜合網絡的前半部分的相關信息,在後半部分它將能夠構建一個圖像

對於那些熟悉ANN應用於機器學習的人來說,這個概念對於自動編碼器的概念非常熟悉,正如Pierre Baldi的 2012年白皮書所解釋的那樣。實際上,綜合相關信息和N到N維度轉換的想法非常相似,但有可能認為這兩種想法都有不同的目標。儘管如此,如果你對自動編碼器很熟悉,那麼很容易就能看到Ronneberger和Christ的論文。

關於ConvNets的更多細節

僅基於卷積和激活函數的純卷積網絡在實際應用中很難收斂到有用的結果。這可能是由於向網絡提供的數據量過大或不足而導致的。

為了解決這些問題,convnet應用諸如池化和dropout等技術。

池化是通過選擇其區域的期望屬性來顯著減少卷積窗口內的數據的方差的操作。例如,可以將卷積區域池化為其平均值,L2範數,距中心像素的加權距離的平均值,或該區域的最大值。池化是一種通常在卷積本身和激活函數之後附加的操作(即:sigmoid,ReLu等)

dropout是另一種旨在避免過度擬合的操作。dropout層以指定的概率隨機地歸零一個整數,因此,至少在理論上,幫助網絡不僅僅依賴於一小部分值,從而避免過度擬合。

機器學習:U-NET ConvNet用於CT掃描分割

目標

這篇文章的目標是使用Ronneberger的U-Net,使用Tensorflow和Keras進行CT掃描,對肝臟和骨骼進行分割。

數據集

為此選擇的CT掃描數據集是3D IRCAD -01(https://www.ircad.fr/research/3dircadb/)。該數據集包含來自20名患者(10名男性和10名女性)的匿名CT掃描。這些圖像由作者以DICOM和VTK格式提供,像素為512x512。對於我們來說,這個數據集的一個缺點是75%的患者出現了肝臟腫瘤,這不是本文的主要分割目標,可能代表了我們試圖分離的數據中的“噪音”。

對於編碼,我們將基本上使用Python,TensorFlow,Keras和Jupyter來應用U-Net。

由於所選數據集中的圖像採用DICOM格式,我們將使用python庫提取此數據,為方便起見,我們將使用SciPy ndimage處理,調整大小和PIL以保存圖像。

DICOM格式

DICOM(醫學數字成像和通信)格式通常用於傳輸醫學圖像數據。除了圖像本身之外,該格式還包含標識和信息標記。DICOM文件中的像素數據可以使用JPEG、JPEG2000、無損JPEG、RLE和LZW格式進行壓縮

由於理解和如何提取圖像數據並不是這個項目的核心,因此選擇了一個Python庫來實現它。

使用的Python庫為dicom,它的安裝很簡單:

pip install dicom

從DICOM文件讀取是通過以下Python代碼實現的:

import dicom

def getDicomFile(filename):

return dicom.read_file(filename)

DICOM格式不僅包含圖像數據,還包含其他信息。

要從DICOM文件獲取圖像數據,可以訪問該pixel_data元素

import dicom

def getPixelDataFromDicom(filename):

return dicom.read_file(filename).pixel_array

3D IRCAD數據集還包含了由醫學專家為這20名患者的所有圖像手工進行的肝臟、骨骼、腫瘤等的真實分割。

機器學習:U-NET ConvNet用於CT掃描分割

上圖左側,來自患者的幾個切片包括在3D IRCARD中。上圖右側,從左側開始在相應切片處對肝臟進行相應的分割。

來自3D IRCAD的DICOM數據是匿名的,但很少有信息(如患者性別和大致年齡)保留下來,官方網頁和下圖所示。對於此特定工作,我們不會使用此信息,但可能與同一數據集中的未來工作相關。

機器學習:U-NET ConvNet用於CT掃描分割

預處理圖像

Ronneberger論文提出的U-Net將572×572像素的圖像作為輸入和輸出。在最低分辨率層,它在flattened之前具有32x32像素。

3D IRCAD的輸入是512x512,但是在個人電腦上訓練這樣的圖像可能需要很長時間才能完成,所以圖像被縮小到128x128像素。

empty pair of images(CT + Mask)由於在訓練過程中不起作用而被丟棄。

我在CT掃描圖像中加入了高斯噪聲,以增加文件的可變性,並處理了黑色背景與圖像感興趣部分的數量關係。值得注意的是,即使是在同一個病人和兩個連續的圖片中應用到圖像上的高斯噪聲也是不同的。

預處理數據集

3D IRCAD,患者數據分佈在兩個文件夾中。在PATIENT_DICOM文件夾中有真實的CT掃描數據,在MASKS_DICOM中有相應的分割文件。

在這兩個文件夾的每一箇中,有20個文件夾編號1到20,代表每個患者。對於MASKS_FOLDER中的每個患者數據,每個分段的生物結構(即:肝臟、骨骼、腫瘤等)

文件夾結構因患者而異,但通常類似於:

機器學習:U-NET ConvNet用於CT掃描分割

在訓練之前,首先我們需要(1)以更簡單的格式轉換圖像以便工作;(2)創建輸入和目標集。

從DICOM格式轉換是通過獲取pixel_array並將其保存為更易處理的格式來完成的。在這種情況下,我選擇了PNG格式。

圖像預處理的Python代碼如下:

import dicom

from scipy.ndimage.interpolation import zoom

from scipy.misc import imsave

def add_gaussian_noise(inp, expected_noise_ratio=0.05):

image = inp.copy()

if len(image.shape) == 2:

row,col= image.shape

ch = 1

else:

row,col,ch= image.shape

mean = 0.

var = 0.1

sigma = var**0.5

gauss = np.random.normal(mean,sigma,(row,col)) * expected_noise_ratio

gauss = gauss.reshape(row,col)

noisy = image + gauss

return noisy

def normalize(img):

arr = img.copy().astype(np.float)

M = np.float(np.max(img))

if M != 0:

arr *= 1./M

return arr

def preprocess(filename, resize_ratio=0.25):

img = dicom.read_file(filename).pixel_array

img = normalize(zoom(img, resize_ratio))

img = add_gaussian_noise(img)

return img

### filelist contains all *.dcm files from 3D IRCAD folder

for dicom_file in filelist:

pp_image = preprocess(dicom_file)

imsave(dicom_file.replace("dcm","png"), pp_image, "png")

通過利用文件夾結構獲得第二預處理任務。在此不描述該方法,因為不需要理解該過程並且它僅基於os.walk方法。

TF + Keras的U-Net

該實驗的最終架構由10個卷積層組成,並且在Dice Score獲得0.87至0.93(取決於我們稍後將探討的一些超參數)。Christ et al. 2016論文使用了19層,Dice Score為0.93 - 0.94。

機器學習:U-NET ConvNet用於CT掃描分割

如上圖所示,為這篇文章實施的最終版本不是2015年Ronneberger等文章的香草版本,也不是2016年Christ等人的版本。主要區別在於每個池化階段後的dropout層和層數。由於來自3D IRCAD的信息量有限,因此選擇了dropout策略。在沒有dropout層的情況下,由於過度擬合,訓練需要更多的交互收斂(如果收斂)。減少層數以加速開發過程,而不是出於技術原因。但是,至少對於上面的數據集和對象,結果非常好。

如上所述,該項目的U-Net由10個卷積層組成。不同的是,從原始作品中添加了一個dropout。

用Keras創建這樣的結構就像,Python代碼如下:

def get_model(optimizer, loss_metric, metrics, lr=1e-3):

inputs = Input((sample_width, sample_height, 1))

conv1 = Conv2D(32, (3, 3), activation='relu', padding='same')(inputs)

conv1 = Conv2D(32, (3, 3), activation='relu', padding='same')(conv1)

pool1 = MaxPooling2D(pool_size=(2, 2))(conv1)

drop1 = Dropout(0.5)(pool1)

conv2 = Conv2D(64, (3, 3), activation='relu', padding='same')(drop1)

conv2 = Conv2D(64, (3, 3), activation='relu', padding='same')(conv2)

pool2 = MaxPooling2D(pool_size=(2, 2))(conv2)

drop2 = Dropout(0.5)(pool2)

conv3 = Conv2D(128, (3, 3), activation='relu', padding='same')(drop2)

conv3 = Conv2D(128, (3, 3), activation='relu', padding='same')(conv3)

pool3 = MaxPooling2D(pool_size=(2, 2))(conv3)

drop3 = Dropout(0.3)(pool3)

conv4 = Conv2D(256, (3, 3), activation='relu', padding='same')(drop3)

conv4 = Conv2D(256, (3, 3), activation='relu', padding='same')(conv4)

pool4 = MaxPooling2D(pool_size=(2, 2))(conv4)

drop4 = Dropout(0.3)(pool4)

conv5 = Conv2D(512, (3, 3), activation='relu', padding='same')(drop4)

conv5 = Conv2D(512, (3, 3), activation='relu', padding='same')(conv5)

up6 = concatenate([Conv2DTranspose(256, (2, 2), strides=(2, 2), padding='same')(conv5), conv4], axis=3)

conv6 = Conv2D(256, (3, 3), activation='relu', padding='same')(up6)

conv6 = Conv2D(256, (3, 3), activation='relu', padding='same')(conv6)

up7 = concatenate([Conv2DTranspose(128, (2, 2), strides=(2, 2), padding='same')(conv6), conv3], axis=3)

conv7 = Conv2D(128, (3, 3), activation='relu', padding='same')(up7)

conv7 = Conv2D(128, (3, 3), activation='relu', padding='same')(conv7)

up8 = concatenate([Conv2DTranspose(64, (2, 2), strides=(2, 2), padding='same')(conv7), conv2], axis=3)

conv8 = Conv2D(64, (3, 3), activation='relu', padding='same')(up8)

conv8 = Conv2D(64, (3, 3), activation='relu', padding='same')(conv8)

up9 = concatenate([Conv2DTranspose(32, (2, 2), strides=(2, 2), padding='same')(conv8), conv1], axis=3)

conv9 = Conv2D(32, (3, 3), activation='relu', padding='same')(up9)

conv9 = Conv2D(32, (3, 3), activation='relu', padding='same')(conv9)

conv10 = Conv2D(1, (1, 1), activation='sigmoid')(conv9)

model = Model(inputs=[inputs], outputs=[conv10])

model.compile(optimizer=optimizer(lr=lr), loss=loss_metric, metrics=metrics)

return model

使用Adam優化器可以獲得最佳結果。Adam Optimizer是Stochastic Gradient Descent的擴展,使用RMSProp和AdaGrad的創意。

Sorosen-Dice係數用於比較兩個不同的統計樣本,但它也用於圖像處理區域以匹配二進制圖像之間的相似性。

損失函數由Dice係數的負數給出。

由於有必要在不同時刻使用Dice係數的計算,考慮到繪圖、記錄和訓練本身,有三種變化。

Python實現如下:

#smooth = 1.

# Dice Coefficient to work with Tensorflow

def dice_coef(y_true, y_pred):

y_true_f = K.flatten(y_true)

y_pred_f = K.flatten(y_pred)

intersection = K.sum(y_true_f * y_pred_f)

return (2. * intersection + smooth) / (K.sum(y_true_f) + K.sum(y_pred_f) + smooth)

def dice_coef_loss(y_true, y_pred):

return -dice_coef(y_true, y_pred)

# Dice Coefficient to work outside Tensorflow

def dice_coef_2(y_true, y_pred):

side = len(y_true[0])

y_true_f = y_true.reshape(side*side)

y_pred_f = y_pred.reshape(side*side)

intersection = sum(y_true_f * y_pred_f)

return (2. * intersection + smooth) / (sum(y_true_f) + sum(y_pred_f) + smooth)

考慮到上面提供的模型和損失函數,我們有訓練模型:

model = get_model(optimizer=Adam, loss_metric=dice_coef_loss, metrics=[dice_coef], lr=1e-3)

訓練epochs的變化和訓練測試的拆分對訓練效果有明顯的影響。但是為了演示,肝臟分割的最佳輸出之一使用了250個epoch,使用了包含感興趣對象的圖片總數的3/4。

機器學習:U-NET ConvNet用於CT掃描分割

結果如上。最左邊的列是CT掃描圖像。中間列是真正的分割,最右邊的列是由我們的10層U-Net生成的分割。

骨分割獲得了類似的結果

機器學習:U-NET ConvNet用於CT掃描分割

結果如上。最左邊的列是CT掃描圖像。中間列是真正的分割,最右邊的列是由我們的10層U-Net生成的分割。

其他結果

我們測試了其他優化器和訓練參數,但沒有更改模型圖。結果顯示如下。

機器學習:U-NET ConvNet用於CT掃描分割

進一步的工作

3D IRCAD提供腫瘤患者的圖像。如果癌症診斷可以完全自動化,則在醫療方面和患者護理成本降低方面取得了相當大的進步。

上面的分割具有相當大的Dice Score,但邊界病例對特定器官檢測(肝臟和骨骼)具有假陽性和假陰性結果。

因為我不是醫生,所以我只能使用數據集提供的ground truth圖像,並假設上面的模型生成了錯誤的信息(對於一個系統來說,最可能的情況是在一個上調整了幾百個已調整大小的圖像或兩兩分鐘)。

增加此分數可能需要使用更多層和使用更大尺寸的圖像來改進模型。這將導致訓練時間成倍增長,這可能因我的個人資源和可用時間而變得不切實際,但在時機成熟時是可能的。

腫瘤的自動化可能不僅僅使用可用的數據。由於腫瘤通常非常大並且通常不具有很多視覺相似性,因此具有20名患者的3D IRCAD數據集可能不足以獲得更好的結果。但僅此一點不應該阻止我們嘗試。


分享到:


相關文章: