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-Nets出现在2015年Ronneberger等人的文章中（https://arxiv.org/abs/1505.04597）。并且在2016年，Christ等人(https://arxiv.org/pdf/1610.02177.pdf)在CT扫描图像上进行自动肝脏分割。关于U-Net的一个很好的想法是，它可以接收图像作为输入，生成另一个图像作为输出，这对于生成分割图像非常方便。
2015年Ronneberger等人展示的U-Net
U-Net的思想是，通过训练，网络将能够通过最小化与所需操作相关的成本函数来综合网络的前半部分的相关信息，在后半部分它将能够构建一个图像
对于那些熟悉ANN应用于机器学习的人来说，这个概念对于自动编码器的概念非常熟悉，正如Pierre Baldi的 2012年白皮书所解释的那样。实际上，综合相关信息和N到N维度转换的想法非常相似，但有可能认为这两种想法都有不同的目标。尽管如此，如果你对自动编码器很熟悉，那么很容易就能看到Ronneberger和Christ的论文。
关于ConvNets的更多细节
仅基于卷积和激活函数的纯卷积网络在实际应用中很难收敛到有用的结果。这可能是由于向网络提供的数据量过大或不足而导致的。
为了解决这些问题，convnet应用诸如池化和dropout等技术。
池化是通过选择其区域的期望属性来显着减少卷积窗口内的数据的方差的操作。例如，可以将卷积区域池化为其平均值，L2范数，距中心像素的加权距离的平均值，或该区域的最大值。池化是一种通常在卷积本身和激活函数之后附加的操作（即：sigmoid，ReLu等）
dropout是另一种旨在避免过度拟合的操作。dropout层以指定的概率随机地归零一个整数，因此，至少在理论上，帮助网络不仅仅依赖于一小部分值，从而避免过度拟合。
目标
这篇文章的目标是使用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名患者的所有图像手工进行的肝脏、骨骼、肿瘤等的真实分割。
上图左侧，来自患者的几个切片包括在3D IRCARD中。上图右侧，从左侧开始在相应切片处对肝脏进行相应的分割。
来自3D IRCAD的DICOM数据是匿名的，但很少有信息(如患者性别和大致年龄)保留下来，官方网页和下图所示。对于此特定工作，我们不会使用此信息，但可能与同一数据集中的未来工作相关。
预处理图像
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中的每个患者数据，每个分段的生物结构(即:肝脏、骨骼、肿瘤等)
文件夹结构因患者而异，但通常类似于:
在训练之前，首先我们需要（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。
如上图所示，为这篇文章实施的最终版本不是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。
结果如上。最左边的列是CT扫描图像。中间列是真正的分割，最右边的列是由我们的10层U-Net生成的分割。
骨分割获得了类似的结果
结果如上。最左边的列是CT扫描图像。中间列是真正的分割，最右边的列是由我们的10层U-Net生成的分割。
其他结果
我们测试了其他优化器和训练参数，但没有更改模型图。结果显示如下。
进一步的工作
3D IRCAD提供肿瘤患者的图像。如果癌症诊断可以完全自动化，则在医疗方面和患者护理成本降低方面取得了相当大的进步。
上面的分割具有相当大的Dice Score，但边界病例对特定器官检测（肝脏和骨骼）具有假阳性和假阴性结果。
因为我不是医生，所以我只能使用数据集提供的ground truth图像，并假设上面的模型生成了错误的信息（对于一个系统来说，最可能的情况是在一个上调整了几百个已调整大小的图像或两两分钟）。
增加此分数可能需要使用更多层和使用更大尺寸的图像来改进模型。这将导致训练时间成倍增长，这可能因我的个人资源和可用时间而变得不切实际，但在时机成熟时是可能的。
肿瘤的自动化可能不仅仅使用可用的数据。由于肿瘤通常非常大并且通常不具有很多视觉相似性，因此具有20名患者的3D IRCAD数据集可能不足以获得更好的结果。但仅此一点不应该阻止我们尝试。
閱讀更多 不靠譜的貓 的文章
