1. 程式人生 > >深度學習(七)U-Net原理以及keras程式碼實現醫學影象眼球血管分割

深度學習(七)U-Net原理以及keras程式碼實現醫學影象眼球血管分割

原文作者:aircraft

原文連結:https://www.cnblogs.com/DOMLX/p/9780786.html

全卷積神經網路

醫學影象分割框架

醫學影象分割主要有兩種框架,一個是基於CNN的,另一個就是基於FCN的。這裡都是通過網路來進行語義分割。

那麼什麼是語義分割?可不是漢字分割句意,在影象處理中有自己的定義。

影象語義分割的意思就是機器自動分割並識別出影象中的內容,比如給出一個人騎摩托車的照片,機器判斷後應當能夠生成右側圖,紅色標註為人,綠色是車(黑色表示 back ground)。

十分鐘看懂影象語義分割技術

所以影象分割對影象理解的意義,就好比讀古書首先要斷句一樣

在 Deeplearning 技術快速發展之前,就已經有了很多做影象分割的技術,其中比較著名的是一種叫做 “Normalized cut” 的圖劃分方法,簡稱 “N-cut

”。

N-cut 的計算有一些連線權重的公式,這裡就不提了,它的思想主要是通過畫素和畫素之間的關係權重來綜合考慮,根據給出的閾值,將影象一分為二

基於CNN 的框架

這個想法也很簡單,就是對影象的每一個畫素點進行分類,在每一個畫素點上取一個patch,當做一幅影象,輸入神經網路進行訓練,舉個例子:

cnnbased

這是一個二分類問題,把影象中所有label為0的點作為負樣本,所有label為1的點作為正樣本。

這種網路顯然有兩個缺點: 1. 冗餘太大,由於每個畫素點都需要取一個patch,那麼相鄰的兩個畫素點的patch相似度是非常高的,這就導致了非常多的冗餘,導致網路訓練很慢。 2. 感受野和定位精度不可兼得,當感受野選取比較大的時候,後面對應的pooling層的降維倍數就會增大,這樣就會導致定位精度降低,但是如果感受野比較小,那麼分類精度就會降低。

基於FCN框架

在醫學影象處理領域,有一個應用很廣泛的網路結構—-U-net ,網路結構如下:

這裡寫圖片描述

可以看出來,就是一個全卷積神經網路,輸入和輸出都是影象,沒有全連線層。較淺的高解析度層用來解決畫素定位的問題,較深的層用來解決畫素分類的問題。

好了理解完U-net網路,我們就學習一下怎麼用U-net網路來進行醫學影象分割。

U-net+kears實現眼部血管分割

原作者的【英文說明】https://github.com/orobix/retina-unet#retina-blood-vessel-segmentation-with-a-convolution-neural-network-u-net

linux下就環境一樣,配置就要自己去找了。

1、介紹為了能夠更好的對眼部血管等進行檢測、分類等操作,我們首先要做的就是對眼底影象中的血管進行分割,保證最大限度的分割出眼部的血管。從而方便後續對血管部分的操作。這部分程式碼選用的資料集是DRIVE資料集,包括訓練集和測試集兩部分。眼底影象資料如圖1所示。

圖1 DRIVE資料集的訓練集眼底影象DRIVE資料集的優點是:不僅有已經手工分好的的血管影象(在manual資料夾下,如圖2所示),而且還包含有眼部輪廓的影象(在mask資料夾下,如圖3所示)。

圖2 DRIVE資料集的訓練集手工標註血管影象

圖3 DRIVE資料集的訓練集眼部輪廓影象DRIVE資料集的缺點是:顯而易見,從上面的圖片中可以看出,訓練集只有20幅圖片,可見資料量實在是少之又少。。。所以,為了得到更好的分割效果,我們需要對這20幅影象進行預處理從而增大其資料量

2、依賴的庫- numpy >= 1.11.1- Keras >= 2.1.0- PIL >=1.1.7- opencv >=2.4.10- h5py >=2.6.0- configparser >=3.5.0b2- scikit-learn >= 0.17.1

3、資料讀取與儲存資料集中訓練集和測試集各只有20幅眼底影象(tif格式)。首先要做的第一步就是對生成資料檔案,方便後續的處理。所以這裡我們需要對資料集中的眼底影象、人工標註的血管影象、眼部輪廓生成資料檔案。這裡使用的是hdf5檔案。有關hdf5檔案的介紹,請參考CSDN部落格(HDF5快速上手全攻略)。

4、網路解析

  因為U-net網路可以針對很少的資料集來進行語義分割,比如我們這個眼球血管分割就是用了20張圖片來訓練就可以達到很好的效果。而且我們這種眼球血管,或者指靜脈,指紋之類的提取特徵或者血管靜脈在U-net網路裡就是一個二分類問題,大家一聽,二分類對於目前的神經網路不是一件很簡單的事情了嗎?還有是什麼可以說的。

  的確目前二分類問題是沒有什麼難度了,只要給我足夠的資料集做訓練。而本文用的U-net網路來實現這個二分類就只需要二十張圖片來作為資料集。大家可以看到優勢所在了吧。

5、具體實現

首先我們肯定都是要對資料進行一些預處理的。

第一步

  先將影象轉為灰度圖分別讀入陣列建立起一個符合我們自己的tensor的格式才好傳入神經網路,這裡我們是先將資料存入hdf5檔案中,在開始執行的時候從檔案中讀入。

#將對應的影象資料存入對應影象陣列
def get_datasets(imgs_dir,groundTruth_dir,borderMasks_dir,train_test="null"):
    imgs = np.empty((Nimgs,height,width,channels))
    groundTruth = np.empty((Nimgs,height,width))
    border_masks = np.empty((Nimgs,height,width))
    for path, subdirs, files in os.walk(imgs_dir): #list all files, directories in the path
        for i in range(len(files)):
            #original
            print("original image: " +files[i])
            img = Image.open(imgs_dir+files[i])
            imgs[i] = np.asarray(img)
            #corresponding ground truth
            groundTruth_name = files[i][0:2] + "_manual1.gif"
            print("ground truth name: " + groundTruth_name)
            g_truth = Image.open(groundTruth_dir + groundTruth_name)
            groundTruth[i] = np.asarray(g_truth)
            #corresponding border masks
            border_masks_name = ""
            if train_test=="train":
                border_masks_name = files[i][0:2] + "_training_mask.gif"
            elif train_test=="test":
                border_masks_name = files[i][0:2] + "_test_mask.gif"
            else:
                print("specify if train or test!!")
                exit()
            print("border masks name: " + border_masks_name)
            b_mask = Image.open(borderMasks_dir + border_masks_name)
            border_masks[i] = np.asarray(b_mask)

    print("imgs max: " +str(np.max(imgs)))
    print("imgs min: " +str(np.min(imgs)))
    assert(np.max(groundTruth)==255 and np.max(border_masks)==255)
    assert(np.min(groundTruth)==0 and np.min(border_masks)==0)
    print("ground truth and border masks are correctly withih pixel value range 0-255 (black-white)")
    #reshaping for my standard tensors
    imgs = np.transpose(imgs,(0,3,1,2))
    assert(imgs.shape == (Nimgs,channels,height,width))
    groundTruth = np.reshape(groundTruth,(Nimgs,1,height,width))
    border_masks = np.reshape(border_masks,(Nimgs,1,height,width))
    assert(groundTruth.shape == (Nimgs,1,height,width))
    assert(border_masks.shape == (Nimgs,1,height,width))
    return imgs, groundTruth, border_masks

if not os.path.exists(dataset_path):
    os.makedirs(dataset_path)
#getting the training datasets
imgs_train, groundTruth_train, border_masks_train = get_datasets(original_imgs_train,groundTruth_imgs_train,borderMasks_imgs_train,"train")
print("saving train datasets")
write_hdf5(imgs_train, dataset_path + "DRIVE_dataset_imgs_train.hdf5")
write_hdf5(groundTruth_train, dataset_path + "DRIVE_dataset_groundTruth_train.hdf5")
write_hdf5(border_masks_train,dataset_path + "DRIVE_dataset_borderMasks_train.hdf5")

#getting the testing datasets
imgs_test, groundTruth_test, border_masks_test = get_datasets(original_imgs_test,groundTruth_imgs_test,borderMasks_imgs_test,"test")
print("saving test datasets")
write_hdf5(imgs_test,dataset_path + "DRIVE_dataset_imgs_test.hdf5")
write_hdf5(groundTruth_test, dataset_path + "DRIVE_dataset_groundTruth_test.hdf5")
write_hdf5(border_masks_test,dataset_path + "DRIVE_dataset_borderMasks_test.hdf5")

第二步

  是對讀入記憶體準備開始訓練的影象資料進行一些增強之類的處理,這裡對其進行了,直方圖均衡化,資料標準化,並且壓縮畫素值到0-1,將其的一個數據符合標準正態分佈。當然啦我們這個資料拿來訓練還是太少的,所以我們對每張圖片取patch時,除了正常的每個patch每個patch移動的取之外,我們還在資料範圍內進行隨機取patch,這樣雖然各個patch之間會有一部分資料是相同的,但是這對於網路而言,你傳入的也是一個新的東西,網路能從中提取到的特徵也更多了。這一步的目的其實就是在有限的資料集中進行一些資料擴充,這也是在神經網路訓練中常用的手段了。

  當然了在這個過程中我們也可以隨機組合小的patch來看看。

隨機原圖:

mask圖:

處理待訓練資料的部分程式碼:

def get_data_training(DRIVE_train_imgs_original,
                      DRIVE_train_groudTruth,
                      patch_height,
                      patch_width,
                      N_subimgs,
                      inside_FOV):
    train_imgs_original = load_hdf5(DRIVE_train_imgs_original)
    train_masks = load_hdf5(DRIVE_train_groudTruth) #masks always the same
    # visualize(group_images(train_imgs_original[0:20,:,:,:],5),'imgs_train')#.show()  #check original imgs train


    train_imgs = my_PreProc(train_imgs_original)    #直方圖均衡化,資料標準化,壓縮畫素值到0-1
    train_masks = train_masks/255.

    train_imgs = train_imgs[:,:,9:574,:]  #cut bottom and top so now it is 565*565
    train_masks = train_masks[:,:,9:574,:]  #cut bottom and top so now it is 565*565
    data_consistency_check(train_imgs,train_masks)

    #check masks are within 0-1
    assert(np.min(train_masks)==0 and np.max(train_masks)==1)

    print("\ntrain images/masks shape:")
    print(train_imgs.shape)
    print("train images range (min-max): " +str(np.min(train_imgs)) +' - '+str(np.max(train_imgs)))
    print("train masks are within 0-1\n")

    #extract the TRAINING patches from the full images
    patches_imgs_train, patches_masks_train = extract_random(train_imgs,train_masks,patch_height,patch_width,N_subimgs,inside_FOV)
    data_consistency_check(patches_imgs_train, patches_masks_train)

    print("\ntrain PATCHES images/masks shape:")
    print(patches_imgs_train.shape)
    print("train PATCHES images range (min-max): " +str(np.min(patches_imgs_train)) +' - '+str(np.max(patches_imgs_train)))

    return patches_imgs_train, patches_masks_train#, patches_imgs_test, patches_masks_test

第三步

按照U-net的網路結構,使用keras來構造出網路。這裡對keras函式語法不太理解的可以看這篇部落格:深度學習(六)keras常用函式學習

def get_unet(n_ch,patch_height,patch_width):
    inputs = Input(shape=(n_ch,patch_height,patch_width))
    conv1 = Conv2D(32, (3, 3), activation='relu', padding='same',data_format='channels_first')(inputs)
    conv1 = Dropout(0.2)(conv1)
    conv1 = Conv2D(32, (3, 3), activation='relu', padding='same',data_format='channels_first')(conv1)
    pool1 = MaxPooling2D((2, 2))(conv1)
    #
    conv2 = Conv2D(64, (3, 3), activation='relu', padding='same',data_format='channels_first')(pool1)
    conv2 = Dropout(0.2)(conv2)
    conv2 = Conv2D(64, (3, 3), activation='relu', padding='same',data_format='channels_first')(conv2)
    pool2 = MaxPooling2D((2, 2))(conv2)
    #
    conv3 = Conv2D(128, (3, 3), activation='relu', padding='same',data_format='channels_first')(pool2)
    conv3 = Dropout(0.2)(conv3)
    conv3 = Conv2D(128, (3, 3), activation='relu', padding='same',data_format='channels_first')(conv3)

    up1 = UpSampling2D(size=(2, 2))(conv3)
    up1 = concatenate([conv2,up1],axis=1)
    conv4 = Conv2D(64, (3, 3), activation='relu', padding='same',data_format='channels_first')(up1)
    conv4 = Dropout(0.2)(conv4)
    conv4 = Conv2D(64, (3, 3), activation='relu', padding='same',data_format='channels_first')(conv4)
    #上取樣後橫向拼接
    up2 = UpSampling2D(size=(2, 2))(conv4)
    up2 = concatenate([conv1,up2], axis=1)
    conv5 = Conv2D(32, (3, 3), activation='relu', padding='same',data_format='channels_first')(up2)
    conv5 = Dropout(0.2)(conv5)
    conv5 = Conv2D(32, (3, 3), activation='relu', padding='same',data_format='channels_first')(conv5)
    #
    conv6 = Conv2D(2, (1, 1), activation='relu',padding='same',data_format='channels_first')(conv5)
    conv6 = core.Reshape((2,patch_height*patch_width))(conv6)
    conv6 = core.Permute((2,1))(conv6)
    ############
    conv7 = core.Activation('softmax')(conv6)

    model = Model(inputs=inputs, outputs=conv7)

    # sgd = SGD(lr=0.01, decay=1e-6, momentum=0.3, nesterov=False)
    model.compile(optimizer='sgd', loss='categorical_crossentropy',metrics=['accuracy'])

    return model

第四步

這裡就是將我們處理好的資料傳入到網路裡訓練了。得出結果圖。

看的出來很多很細的紋理都被提取出來了,這個U-net網路也可以用於一些醫學細胞的 邊緣提取,指靜脈,掌靜脈之類的紋路提取都可以。後面在下可能還會出指靜脈之類其他影象的語義分割提取,關注在下就可以看到啦hhhhhhh

因為keras內部可以直接將整個網路結構打印出來,所以我們可以看到完整的網路結構圖如下所示:

參考部落格:https://blog.csdn.net/u013063099/article/details/79981097

參考部落格:https://blog.csdn.net/qq_16900751/article/details/78251778