在之前的一篇博客里:医学图像分割 unet实现(一),是学习并复现别人的实验。这篇将记录下自己毕设第一阶段的实验。毕设题目为:基于深度学习的肝脏肿瘤分割。
经过几番调整,最终确定:第一阶段分割出腹部图像中的肝脏,作为第二阶段的ROI(region of interest),第二阶段利用ROI对腹部图像进行裁剪,裁剪后的非ROI区域变成黑色,作为该阶段输入,分割出肝脏中的肿瘤(更新2019-4-2,已做完实验,医学图像分割 基于深度学习的肝脏肿瘤分割 实战(二)。第三阶段用随机场的后处理方法进行优化。
此外,觉得当时的那篇文章写得很不用心,没费多大功夫,复制粘贴为主。我记录的主要原因之一就是觉得这方面入门实战的文章很少,希望让情况类似的同学少走一点弯路。可现在自己写的东西过几天再去看,都会一头雾水。没有从主干思路到细节的分层讲解,没有关键点的介绍。这篇文章开始,将尽我所能,写得利于接受一点。
(不过,自己才开始学习一个多月,python和框架都是临时学的,并不能保证博客里的东西是对的,有错误的话,望指出)

正文:

文章目录

  • 目标
  • 原始数据介绍
  • 整体思路
    • 1、数据提取
    • 2、数据的预处理
    • 3、数据增强
    • 4、数据存储
    • 测试数据存储的全部过程
    • 5、构建网络
    • 6、进行训练并测试

目标

分割出CT腹部图像的肝脏区域。

原始数据介绍

实验用到的数据为3Dircadb,即腹部的CT图。一个病人为一个文件夹。
例如3Dircadb1.1(一共20位),该文件夹下会使用到的子文件夹为PATIENT_DICOM(病人原始腹部3D CT图),MASKS_DICOM(该文件夹下具有不同部位的分割结果,比如liver和liver tumor等等)。如下图所示:

PATIENT_DICOM利用软件展示效果如下:一个dcm文件包含129张切片。

MASKS_DICOM下的liver分割图效果如下:

整体思路

1、数据提取

数据读取:
从原始dcm格式读入成我们需要的数组格式

#part1
import numpy as np
import pydicom
import os
import matplotlib.pyplot as plt
import cv2
from keras.preprocessing.image import ImageDataGenerator
from HDF5DatasetWriter import HDF5DatasetWriter
from HDF5DatasetGenerator import HDF5DatasetGeneratorfor i in range(1,18): # 前17个人作为测试集full_images = [] # 后面用来存储目标切片的列表full_livers = [] #功能同上# 注意不同的系统,文件分割符的区别label_path = '~/3Dircadb/3Dircadb1.%d/MASKS_DICOM/liver'%idata_path = '~/3Dircadb/3Dircadb1.%d/PATIENT_DICOM'%iliver_slices = [pydicom.dcmread(label_path + '/' + s) for s in os.listdir(label_path)]# 注意需要排序,即使文件夹中显示的是有序的,读进来后就是随机的了liver_slices.sort(key = lambda x: int(x.InstanceNumber))# s.pixel_array 获取dicom格式中的像素值livers = np.stack([s.pixel_array for s in liver_slices])image_slices = [pydicom.dcmread(data_path + '/' + s) for s in os.listdir(data_path)]image_slices.sort(key = lambda x: int(x.InstanceNumber))""" 省略进行的预处理操作,具体见part2"""full_images.append(images)full_livers.append(livers)full_images = np.vstack(full_images)full_images = np.expand_dims(full_images,axis=-1)full_livers = np.vstack(full_livers)full_livers = np.expand_dims(full_livers,axis=-1)

2、数据的预处理

1、将ct值转化为标准的hu值
至于为什么要将值进行转化,这儿就不详细说明,具体可以参考医学图像预处理(三)——windowing(ct对比增强),这篇博文中有一样的get_pixels_hu函数
2、窗口化操作
医学图像预处理(三)——windowing(ct对比增强)
3、直方图均衡化

def clahe_equalized(imgs,start,end):assert (len(imgs.shape)==3)  #3D arrays#create a CLAHE object (Arguments are optional).clahe = cv2.createCLAHE(clipLimit=2.0, tileGridSize=(8,8))imgs_equalized = np.empty(imgs.shape)for i in range(start, end+1):imgs_equalized[i,:,:] = clahe.apply(np.array(imgs[i,:,:], dtype = np.uint8))return imgs_equalized

4、归一化
5、仅提取腹部所有切片中包含了肝脏的那些切片,其余的不要
医学图像预处理(四)—— 提取包含目标的切片

#part2# 接part1images = get_pixels_hu(image_slices)images = transform_ctdata(images,500,150)start,end = getRangImageDepth(livers)images = clahe_equalized(images,start,end)images /= 255.# 仅提取腹部所有切片中包含了肝脏的那些切片,其余的不要total = (end - 4) - (start+4) +1print("%d person, total slices %d"%(i,total))# 首和尾目标区域都太小,舍弃images = images[start+5:end-5]print("%d person, images.shape:(%d,)"%(i,images.shape[0]))livers[livers>0] = 1livers = livers[start+5:end-5]

3、数据增强

利用keras的数据增强接口,可以实现分割问题的数据增强。一般的增强是分类问题,这种情况,只需要对image变形,label保持不变。但分割问题,就需要image和mask进行同样的变形处理。具体怎么实现,参考下面代码,注意种子设定成一样的。

# 可以在part1之前设定好(即循环外)
seed=1
data_gen_args = dict(rotation_range=3,width_shift_range=0.01,height_shift_range=0.01,shear_range=0.01,zoom_range=0.01,fill_mode='nearest')image_datagen = ImageDataGenerator(**data_gen_args)
mask_datagen = ImageDataGenerator(**data_gen_args)
#part3 接part2image_datagen.fit(full_images, augment=True, seed=seed)mask_datagen.fit(full_livers, augment=True, seed=seed)image_generator = image_datagen.flow(full_images,seed=seed)mask_generator = mask_datagen.flow(full_livers,seed=seed)train_generator = zip(image_generator, mask_generator)x=[]y=[]i = 0for x_batch, y_batch in train_generator:i += 1x.append(x_batch)y.append(y_batch)if i>=2: # 因为我不需要太多的数据breakx = np.vstack(x)y = np.vstack(y)

4、数据存储

一般而言,数据量较大的话,都会先将原始数据库的东西转化为np或者h5格式的文件,我感觉这样有两个好处,一是真正输入网络训练的时候io量会大大减少(特别是h5很适用于大的数据库),二是数据分享或者上传至服务器时也方便一点。

实验中会出现两个类,分别是写h5和读h5文件的辅助类
这读文件的类写成了generator,这样可以结合训练网络时,keras的fit_generator来使用,降低内存开销。

class HDF5DatasetWriter:"""用来写数据到h5文件"""
class HDF5DatasetGenerator:"""用来读h5文件的数据"""

它们具体的实现在python h5文件的读写,因为篇幅问题,所以这儿不详述

h5文件操作需要的包import h5py

# 可以在part1之前设定好(即循环外)
# 这儿的数量需要提前写好,感觉很不方便,但不知道怎么改,我是先跑了之前的程序,计算了一共有多少
# 张图片后再写的,但这样明显不是好的解决方案
dataset = HDF5DatasetWriter(image_dims=(2782, 512, 512, 1),mask_dims=(2782, 512, 512, 1),outputPath="data_train/train_liver.h5")
#part4 接part3dataset.add(full_images, full_livers)dataset.add(x, y)# end of lop
dataset.close()

测试数据存储的全部过程

测试数据与上面训练数据的处理过程几乎一样,但测试数据不要进行数据增强

full_images2 = []
full_livers2 = []
for i in range(18,21):#后3个人作为测试样本label_path = '~/3Dircadb/3Dircadb1.%d/MASKS_DICOM/liver'%idata_path = '~/3Dircadb/3Dircadb1.%d/PATIENT_DICOM'%iliver_slices = [pydicom.dcmread(label_path + '/' + s) for s in os.listdir(label_path)]liver_slices.sort(key = lambda x: int(x.InstanceNumber))livers = np.stack([s.pixel_array for s in liver_slices])start,end = getRangImageDepth(livers)total = (end - 4) - (start+4) +1print("%d person, total slices %d"%(i,total))image_slices = [pydicom.dcmread(data_path + '/' + s) for s in os.listdir(data_path)]image_slices.sort(key = lambda x: int(x.InstanceNumber))images = get_pixels_hu(image_slices)images = transform_ctdata(images,500,150)images = clahe_equalized(images,start,end)images /= 255.images = images[start+5:end-5]print("%d person, images.shape:(%d,)"%(i,images.shape[0]))livers[livers>0] = 1livers = livers[start+5:end-5]full_images2.append(images)full_livers2.append(livers)full_images2 = np.vstack(full_images2)
full_images2 = np.expand_dims(full_images2,axis=-1)
full_livers2 = np.vstack(full_livers2)
full_livers2 = np.expand_dims(full_livers2,axis=-1)dataset = HDF5DatasetWriter(image_dims=(full_images2.shape[0], full_images2.shape[1], full_images2.shape[2], 1),mask_dims=(full_images2.shape[0], full_images2.shape[1], full_images2.shape[2], 1),outputPath="data_train/val_liver.h5")dataset.add(full_images2, full_livers2)print("total images in val ",dataset.close())

5、构建网络

这一部分不多说,是直接用的别人写好的Unet。
唯一进行的更改,就是Crop的那部分。因为,如果输入图片的高或者宽不能等于2^n(n>=m),m为网络收缩路径的层数。那么出现的问题就是,下采样的时候进行了取整操作,但是,后续在扩张路径又会有上采样,且上采样的结果会和收缩路径的特征图进行拼接,即long skip connection。若没有达到上面提到的要求,会导致拼接的两个特征图size不同,报错。

# partA
import os
import sys
import numpy as np
import random
import math
import tensorflow as tf
from HDF5DatasetGenerator import HDF5DatasetGenerator
from keras.models import Model
from keras.layers import Input, concatenate, Conv2D, MaxPooling2D, Conv2DTranspose,Cropping2D
from keras.optimizers import Adam
from keras.callbacks import ModelCheckpoint
from keras import backend as K
from skimage import ioK.set_image_data_format('channels_last')def dice_coef(y_true, y_pred):smooth = 1.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)def get_crop_shape(target, refer):# width, the 3rd dimensionprint(target.shape)print(refer._keras_shape)cw = (target._keras_shape[2] - refer._keras_shape[2])assert (cw >= 0)if cw % 2 != 0:cw1, cw2 = int(cw/2), int(cw/2) + 1else:cw1, cw2 = int(cw/2), int(cw/2)# height, the 2nd dimensionch = (target._keras_shape[1] - refer._keras_shape[1])assert (ch >= 0)if ch % 2 != 0:ch1, ch2 = int(ch/2), int(ch/2) + 1else:ch1, ch2 = int(ch/2), int(ch/2)return (ch1, ch2), (cw1, cw2)def get_unet():inputs = Input((IMG_HEIGHT, IMG_WIDTH , 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)conv2 = Conv2D(64, (3, 3), activation='relu', padding='same')(pool1)conv2 = Conv2D(64, (3, 3), activation='relu', padding='same')(conv2)pool2 = MaxPooling2D(pool_size=(2, 2))(conv2)conv3 = Conv2D(128, (3, 3), activation='relu', padding='same')(pool2)conv3 = Conv2D(128, (3, 3), activation='relu', padding='same')(conv3)pool3 = MaxPooling2D(pool_size=(2, 2))(conv3)conv4 = Conv2D(256, (3, 3), activation='relu', padding='same')(pool3)conv4 = Conv2D(256, (3, 3), activation='relu', padding='same')(conv4)pool4 = MaxPooling2D(pool_size=(2, 2))(conv4)conv5 = Conv2D(512, (3, 3), activation='relu', padding='same')(pool4)conv5 = Conv2D(512, (3, 3), activation='relu', padding='same')(conv5)up_conv5 = Conv2DTranspose(256, (2, 2), strides=(2, 2), padding='same')(conv5)ch, cw = get_crop_shape(conv4, up_conv5)crop_conv4 = Cropping2D(cropping=(ch,cw), data_format="channels_last")(conv4)up6 = concatenate([up_conv5, crop_conv4], axis=3)conv6 = Conv2D(256, (3, 3), activation='relu', padding='same')(up6)conv6 = Conv2D(256, (3, 3), activation='relu', padding='same')(conv6)up_conv6 = Conv2DTranspose(128, (2, 2), strides=(2, 2), padding='same')(conv6)ch, cw = get_crop_shape(conv3, up_conv6)crop_conv3 = Cropping2D(cropping=(ch,cw), data_format="channels_last")(conv3)up7 = concatenate([up_conv6, crop_conv3], axis=3)conv7 = Conv2D(128, (3, 3), activation='relu', padding='same')(up7)conv7 = Conv2D(128, (3, 3), activation='relu', padding='same')(conv7)up_conv7 = Conv2DTranspose(64, (2, 2), strides=(2, 2), padding='same')(conv7)ch, cw = get_crop_shape(conv2, up_conv7)crop_conv2 = Cropping2D(cropping=(ch,cw), data_format="channels_last")(conv2)up8 = concatenate([up_conv7, crop_conv2], axis=3)conv8 = Conv2D(64, (3, 3), activation='relu', padding='same')(up8)conv8 = Conv2D(64, (3, 3), activation='relu', padding='same')(conv8)up_conv8 = Conv2DTranspose(64, (2, 2), strides=(2, 2), padding='same')(conv8)ch, cw = get_crop_shape(conv1, up_conv8)crop_conv1 = Cropping2D(cropping=(ch,cw), data_format="channels_last")(conv1)up9 = concatenate([up_conv8, crop_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=Adam(lr=1e-5), loss=dice_coef_loss, metrics=[dice_coef])return model

6、进行训练并测试

这其中主要包括训练文件与测试文件的读取,Checkpoint回掉函数的设定,fit_generator的使用,模型预测,并对预测结果进行保存。

# partB 接partA
IMG_WIDTH = 512
IMG_HEIGHT = 512
IMG_CHANNELS = 1
TOTAL = 2782 # 总共的训练数据
TOTAL_VAL = 152 # 总共的validation数据
# part1部分储存的数据文件
outputPath = './data_train/train_liver.h5' # 训练文件
val_outputPath = './data_train/val_liver.h5'
#checkpoint_path = 'model.ckpt'
BATCH_SIZE = 8 # 根据服务器的GPU显存进行调整class UnetModel:def train_and_predict(self):reader = HDF5DatasetGenerator(dbPath=outputPath,batchSize=BATCH_SIZE)train_iter = reader.generator()test_reader = HDF5DatasetGenerator(dbPath=val_outputPath,batchSize=BATCH_SIZE)test_iter = test_reader.generator()fixed_test_images, fixed_test_masks = test_iter.__next__()
#   model = get_unet()model_checkpoint = ModelCheckpoint('weights.h5', monitor='val_loss', save_best_only=True)# 注:感觉validation的方式写的不对,应该不是这样弄的model.fit_generator(train_iter,steps_per_epoch=int(TOTAL/BATCH_SIZE),verbose=1,epochs=500,shuffle=True,validation_data=(fixed_test_images, fixed_test_masks),callbacks=[model_checkpoint])
#        reader.close()test_reader.close()print('-'*30)print('Loading and preprocessing test data...')print('-'*30)print('-'*30)print('Loading saved weights...')print('-'*30)model.load_weights('weights.h5')print('-'*30)print('Predicting masks on test data...')print('-'*30)imgs_mask_test = model.predict(fixed_test_images, verbose=1)np.save('imgs_mask_test.npy', imgs_mask_test)print('-' * 30)print('Saving predicted masks to files...')print('-' * 30)pred_dir = 'preds'if not os.path.exists(pred_dir):os.mkdir(pred_dir)i = 0for image in imgs_mask_test:image = (image[:, :, 0] * 255.).astype(np.uint8)gt = (fixed_test_masks[i,:,:,0] * 255.).astype(np.uint8)ini = (fixed_test_images[i,:,:,0] *255.).astype(np.uint8)io.imsave(os.path.join(pred_dir, str(i) + '_ini.png'), ini)io.imsave(os.path.join(pred_dir, str(i) + '_pred.png'), image)io.imsave(os.path.join(pred_dir, str(i) + '_gt.png'), gt)i += 1unet = UnetModel()
unet.train_and_predict()

模型跑的过程如图。

预测结果可视化展示

医学图像分割 基于深度学习的肝脏肿瘤分割 实战(一)相关推荐

  1. 医学图像分割 基于深度学习的肝脏肿瘤分割 实战(二)

    在医学图像分割 基于深度学习的肝脏肿瘤分割 实战(一)中,实现了对肝脏的分割,但是后续在使用相同的处理方法与模型进行肿瘤分割的时候,遇到了两次问题. 第一次,网络的dice系数,训练集上一直只能达到4 ...

  2. (脑肿瘤分割笔记:四十四)基于深度学习的脑肿瘤分割的综述

    目录 Abstract&Introduction 脑肿瘤分割任务面临的主要挑战 深度学习方法的脑肿瘤分割的方法 脑肿瘤分割方法一:设计有效的架构分割方法 针对于不同目的的模型 对于精度有要求的 ...

  3. 基于深度学习的点云分割网络及点云分割数据集

    作者丨泡椒味的泡泡糖 来源丨深蓝AI 引言 点云分割是根据空间.几何和纹理等特征对点云进行划分,使得同一划分内的点云拥有相似的特征.点云的有效分割是许多应用的前提,例如在三维重建领域,需要对场景内的物 ...

  4. 前沿丨基于深度学习的点云分割网络及点云分割数据集

    众所周知,点云的有效分割是许多应用的前提,例如在三维重建领域,需要对场景内的物体首先进行分类处理,然后才能进行后期的识别和重建.传统的点云分割主要依赖聚类算法和基于随机采样一致性的分割算法,在很多技术 ...

  5. 大盘点 | 基于深度学习的点云分割网络及点云分割数据集

    编辑 | 深蓝前沿 点击下方卡片,关注"自动驾驶之心"公众号 ADAS巨卷干货,即可获取 后台回复[数据集下载]获取计算机视觉近30种数据集! 引言 点云分割是根据空间.几何和纹理 ...

  6. 语义分割源代码_综述 | 基于深度学习的实时语义分割方法:全面调研

    34页综述,共计119篇参考文献.本文对图像分割中的最新深度学习体系结构进行了全面分析,更重要的是,它提供了广泛的技术列表以实现快速推理和计算效率. A Survey on Deep Learning ...

  7. 基于深度学习的图像语义分割技术概述之4常用方法

    本文为论文阅读笔记,不当之处,敬请指正. A Review on Deep Learning Techniques Applied to Semantic Segmentation:原文链接 4 深度 ...

  8. 基于深度学习的图像语义分割技术概述之背景与深度网络架构

    本文为论文阅读笔记,不当之处,敬请指正.  A Review on Deep Learning Techniques Applied to Semantic Segmentation: 原文链接 摘要 ...

  9. 医学图像分割的深度学习:综述

    Title:Deep Learning for Medical Image Segmentation:Tricks, Challenges and Future Directions 摘要 最近的医学 ...

最新文章

  1. iOS完全自学手册——[三]Objective-C语言速成,利用Objective-C创建自己的对象...
  2. Java语言 泛型 类型擦除
  3. 【干货】功能堆砌or视觉美观?看优秀PM如何权衡
  4. 鸡蛋学运维-2:Rsync同步配置步骤
  5. mysql 超长记录_谁记录了mysql error log中的超长信息(记pt-stalk一个bug的定位过程)...
  6. fedora 不在sudoers文件中_COPR 仓库中 4 个很酷的新软件(2019.4) | Linux 中国
  7. 南京林业大学计算机专业分数线,2021南京林业大学录取分数线_历年各专业分数线(2017-2020),各省投档线_一品高考网...
  8. 大数据_Hbase-内容回顾和补充---Hbase工作笔记0018
  9. (cljs/run-at (JSVM. :browser) 简单类型可不简单啊~)
  10. 01Hypertext Preprocessor
  11. 输出重定向与输入重定向
  12. Linux权限与sudo
  13. Android icon图标网站
  14. GOP和帧率、码率的关系
  15. 海归35岁,阿里P7offer, 是否接受?
  16. windows xp快捷键
  17. 冲击港交所:百果园书写水果连锁运营默示录
  18. C#计算wgs84大地坐标转换为空间直角坐标
  19. 基于 MQTT 通讯一个简单的 Java工程
  20. python mysql_config not found_解决问题:OSError: mysql_config not found

热门文章

  1. 各种幻灯片切换效果。soChange
  2. kaggle注册验证码问题
  3. 文本相似度之五种无监督算法实现代码
  4. 台式机ubuntu,使用intel核显作显示输出,nvidia独显做cuda运算
  5. 安全、开放、智能,华为云WeLink应对云上办公新挑战
  6. 通信屌丝也谈星际穿越,通信大牛请绕道
  7. 【创新基金项目可行性报告】如何写项目可行性报告/项目申报文档
  8. Cocos2d-js中Chipmunk引擎
  9. 无纸化自动化考评系统,专为企业研发安全稳定
  10. 关于ADRC应用在四旋翼飞行器上的论文还挺多的