DLTK是用于医学成像的深度学习工具包,它扩展了TensorFlow, 以实现生物医学图像的深度学习。它提供了专业操作和功能,模型实现,教程(在本博客中使用)和典型应用程序的代码示例。
这篇文章是对生物医学图像深度学习的快速介绍,我们将在这里展示当前工程问题的一些问题和解决方案,并向您展示如何启动和运行问题的原型。
概述
什么是生物医学图像分析?为什么需要它?
生物医学图像是在不同尺度上(即微观,宏观等)对人体的测量。它们具有多种成像模式(例如CT扫描仪,超声波仪器等)并且测量人体的物理特性(例如,辐射密度,X射线的不透明度)。这些图像由领域专家(例如放射科医师)解释用于临床任务(例如诊断)并且对医生的决策具有很大影响。
医学图像示例
生物医学图像通常是体积图像(3D),并且有时具有额外的时间维度(4D)和/或多个通道(4-5D)(例如,多序列MR图像)。生物医学图像的变化与自然图像(例如照片)的变化完全不同,因为临床方案旨在对图像的获取方式进行分层(如病人仰卧,头部不倾斜等)。在他们的分析中,我们的目标是检测细微的差异(即一些小区域表明发现异常)。
为什么需要计算机视觉和机器学习?
长期以来,计算机视觉方法被用于自动分析生物医学图像。最近深度学习的出现取代了许多其他机器学习方法,因为它避免了手工工程特征的创建,从而从过程中消除了一个关键的错误来源。此外,GPU加速的完全网络的快速推理速度允许我们将分析量化为前所未有的数据量(例如10⁶主题图片)。
我们可以随时使用深度学习库进行生物医学成像吗?为什么要创建DLTK?
创建DLTK的主要原因是将这个领域的专用工具开箱即用。虽然许多深度学习库公开低级操作(如张量乘法等)开发人员,但许多较高级别的专业操作因其在体积图像上的使用而丢失(如differentiable 3D upsampling layers层,),并且由于图像的额外空间维度,我们可能遇到存储器问题(例如,存储1k CT图像的数据库的单个副本,其中图像尺寸为325x512x256在float32中为~268GB)。由于采集的性质不同,由于采集的性质不同一些图像将需要特殊的预处理(如强度归一化、偏差校正、去噪、空间归一化/registration等)。
文件格式,headers和读取图像
虽然许多成像模式的供应商都以DICOM标准格式生成图像,以二维切片的形式保存大量数据,但许多分析库都依赖于更适合计算和与医学图像交互的格式。我们使用最初为脑成像开发的NifTI(或.nii格式),但在DLTK和本教程中广泛应用于大多数其他卷图像。这种和其他格式保存的是重建图像容器并将其定位在物理空间的必要信息。
为此,它需要专业的标题信息,我们将通过一些属性进行深入学习:
- 尺寸和大小存储有关如何重建图像的信息(例如,使用大小向量将体积分解为三维)。
- 数据类型
- 体素间距(也是体素的物理尺寸,通常以mm为单位)
- 物理坐标系原点
- 方向
为什么这些属性很重要?
网络将在体素空间中进行训练,这意味着我们将创建形状和尺寸张量[batch_size,dx,dy,dz,channels / features]并将其提供给网络。网络将在该体素空间中进行训练并假设所有图像(也是看不见的测试图像)在该空间中被标准化或者可能有问题要概括。在该体素空间中,特征提取器(例如卷积层)将假设体素尺寸是各向同性的(即,在每个维度中相同)并且所有图像以相同的方式定向。
但是,由于大多数图像都描绘了物理空间,我们需要从该物理空间转换为常见的体素空间:
如果所有图像都以相同的方式定向(有时我们需要注册以对图像进行空间标准化:检查MIRTK),我们可以计算从物理空间到体素空间的缩放变换
phys_coords = origin voxel_spacing * voxel_coord
其中所有这些信息都是存储在.nii头文件中的向量。
读取.nii图像:
有几个库可以读取.nii文件并访问头信息并解析它以获得重建的图像容器作为numpy数组。我们选择了SimpleITK,它允许我们导入额外的图像过滤器以进行预处理和其他任务:
import SimpleITK as sitk
import numpy as np
# A path to a T1-weighted brain .nii image:
t1_fn = './brain_t1_0001.nii'
# Read the .nii image containing the volume with SimpleITK:
sitk_t1 = sitk.ReadImage(t1_fn)
# and access the numpy array:
t1 = sitk.GetArrayFromImage(sitk_t1)
数据I / O方面的考虑
根据训练数据库的大小,有几个选项可以将.nii图像数据输入到网络图中。每一种方法在速度上都有特定的权衡,在训练过程中都可能成为瓶颈。我们将讨论并解释三种选择:
1、In memory & feeding dictionaries
我们可以为网络图创建一个tf.placeholder,并在训练期间通过feed_dict提供它。我们从磁盘读取所有.nii文件,在python中处理它们(cf load_data())并将所有训练示例存储在内存中:
# Load all data into memory
data = load_data(all_filenames, tf.estimator.ModeKeys.TRAIN, reader_params)
# Create placeholder variables and define their shapes (here,
# we input a volume image of size [128, 224, 244] and a single
# channel (i.e. greyscale):
x = tf.placeholder(reader_example_dtypes['features']['x'],
[None, 128, 224, 224, 1])
y = tf.placeholder(reader_example_dtypes['labels']['y'],
[None, 1])
# Create a tf.data.Dataset
dataset = tf.data.Dataset.from_tensor_slices((x, y))
dataset = dataset.repeat(None)
dataset = dataset.batch(batch_size)
dataset = dataset.prefetch(1)
# Create an iterator
iterator = dataset.make_initializable_iterator()
nx = iterator.get_next()
with tf.train.MonitoredTrainingSession() as sess_dict:
sess_dict.run(iterator.initializer,
feed_dict={x: data['features'], y: data['labels']})
for i in range(iterations):
# Get next features/labels pair
dict_batch_feat, dict_batch_lbl = sess_dict.run(nx)
TLDR:这种直接的方法通常是最快且最容易实现的,因为它可以避免连续读取磁盘上的数据,但是需要在内存中保留整个训练示例(和验证示例)的数据库,这对于大型数据库或较大的图像文件来说是不可行的。
2、使用TFRecords数据库:
对于图像卷上的大多数深度学习问题,训练示例的数据库太大,无法装入内存。TFRecords格式允许序列化训练实例并将其存储在磁盘上,用快速写访问(即并行数据读取)
def _int64_feature(value):
return tf.train.Feature(int64_list=tf.train.Int64List(value=[value]))
def _float_feature(value):
return tf.train.Feature(float_list=tf.train.FloatList(value=value))
# path to save the TFRecords file
train_filename = 'train.tfrecords'
# open the file
writer = tf.python_io.TFRecordWriter(train_filename)
# iterate through all .nii files:
for meta_data in all_filenames:
# Load the image and label
img, label = load_img(meta_data, reader_params)
# Create a feature
feature = {'train/label': _int64_feature(label),
'train/image': _float_feature(img.ravel())}
# Create an example protocol buffer
example = tf.train.Example(features=tf.train.Features(feature=feature))
# Serialize to string and write on the file
writer.write(example.SerializeToString())
writer.close()
格式可以直接与TensorFlow连接,可以直接集成到tf.graph中的训练循环中
def decode(serialized_example):
# Decode examples stored in TFRecord
# NOTE: make sure to specify the correct dimensions for the images
features = tf.parse_single_example(
serialized_example,
features={'train/image': tf.FixedLenFeature([128, 224, 224, 1], tf.float32),
'train/label': tf.FixedLenFeature([], tf.int64)})
# NOTE: No need to cast these features, as they are already `tf.float32` values.
return features['train/image'], features['train/label']
dataset = tf.data.TFRecordDataset(train_filename).map(decode)
dataset = dataset.repeat(None)
dataset = dataset.batch(batch_size)
dataset = dataset.prefetch(1)
iterator = dataset.make_initializable_iterator()
features, labels = iterator.get_next()
nx = iterator.get_next()
with tf.train.MonitoredTrainingSession() as sess_rec:
sess_rec.run(iterator.initializer)
for i in range(iterations):
try:
# Get next features-labels pair
rec_batch_feat, rec_batch_lbl = sess_rec.run([features, labels])
except tf.errors.OutOfRangeError:
pass
TLDR : TFRecords是从磁盘访问文件的快速方法,但需要存储整个训练数据库的另一个副本。如果我们的目标是使用几TB大小的数据库,这可能会让人望而却步。
3、使用本机python生成器
最后,我们可以使用python生成器,创建一个read_fn()来直接加载图像数据......
def read_fn(file_references, mode, params=None):
# We define a `read_fn` and iterate through the `file_references`, which
# can contain information about the data to be read (e.g. a file path):
for meta_data in file_references:
# Here, we parse the `subject_id` to construct a file path to read
# an image from.
subject_id = meta_data[0]
data_path = '../../data/IXI_HH/1mm'
t1_fn = os.path.join(data_path, '{}/T1_1mm.nii.gz'.format(subject_id))
# Read the .nii image containing a brain volume with SimpleITK and get
# the numpy array:
sitk_t1 = sitk.ReadImage(t1_fn)
t1 = sitk.GetArrayFromImage(sitk_t1)
# Normalise the image to zero mean/unit std dev:
t1 = whitening(t1)
# Create a 4D Tensor with a dummy dimension for channels
t1 = t1[..., np.newaxis]
# If in PREDICT mode, yield the image (because there will be no label
# present). Additionally, yield the sitk.Image pointer (including all
# the header information) and some metadata (e.g. the subject id),
# to facilitate post-processing (e.g. reslicing) and saving.
# This can be useful when you want to use the same read function as
# python generator for deployment.
if mode == tf.estimator.ModeKeys.PREDICT:
yield {'features': {'x': t1}}
# Labels: Here, we parse the class *sex* from the file_references
# \in [1,2] and shift them to \in [0,1] for training:
sex = np.int32(meta_data[1]) - 1
y = sex
# If training should be done on image patches for improved mixing,
# memory limitations or class balancing, call a patch extractor
if params['extract_examples']:
images = extract_random_example_array(
t1,
example_size=params['example_size'],
n_examples=params['n_examples'])
# Loop the extracted image patches and yield
for e in range(params['n_examples']):
yield {'features': {'x': images[e].astype(np.float32)},
'labels': {'y': y.astype(np.int32)}}
# If desired (i.e. for evaluation, etc.), return the full images
else:
yield {'features': {'x': images},
'labels': {'y': y.astype(np.int32)}}
return
tf.data.Dataset.from_generator()对示例进行排队:
# Generator function
def f():
fn = read_fn(file_references=all_filenames,
mode=tf.estimator.ModeKeys.TRAIN,
params=reader_params)
ex = next(fn)
# Yield the next image
yield ex
# Timed example with generator io
dataset = tf.data.Dataset.from_generator(
f, reader_example_dtypes, reader_example_shapes)
dataset = dataset.repeat(None)
dataset = dataset.batch(batch_size)
dataset = dataset.prefetch(1)
iterator = dataset.make_initializable_iterator()
next_dict = iterator.get_next()
with tf.train.MonitoredTrainingSession() as sess_gen:
# Initialize generator
sess_gen.run(iterator.initializer)
with Timer('Generator'):
for i in range(iterations):
# Fetch the next batch of images
gen_batch_feat, gen_batch_lbl = sess_gen.run([next_dict['features'], next_dict['labels']])
TLDR :它避免创建图像数据库的其他副本,但是比TFRecords慢得多,因为生成器无法并行读取和映射函数。
速度基准和选择方法:我们将这三种读取.nii文件到TensorFlow的方法进行了运行,并比较了加载和提供固定大小的示例数据库所需的时间。
显然,最快的方法是在5.6秒内通过占位符从内存中获取信息,然后是31.1秒的TFRecords,以及123.5秒的未优化磁盘读取。然而,只要训练过程中的forward/backward passes是计算瓶颈,数据I/O的速度就可以忽略不计。
数据归一化
与自然图像一样,我们可以对生物医学图像数据进行归一化,但是方法可能略有不同。归一化的目的是消除数据中的某些变化(例如,不同的主题或图像对比度等),这是已知的,因此简化了我们感兴趣的细微差异的检测(例如,病理学的存在)。在这里,我们将讨论最常见的归一化形式:
体素强度归一化:这种形式高度依赖于成像方式,数据是已获得的。典型的零均值、单位方差归一化是定性图像(例如加权脑MR图像,其对比度高度依赖于采集参数,通常由专家设置)的标准。如果我们采用这种统计方法,我们使用的统计数据来自一个完整的卷,而不是整个数据库。
与此相反,定量成像测量的是物理量(例如CT成像中的辐射密度,在不同的扫描仪中强度相当),并受益于剪切和/或重新缩放,如简单范围归一化(例如到[-1,1])。
实例:强度归一化方法
空间归一化:对图像方向进行归一化避免了模型必须学习所有可能的方向,这大大减少了所需的训练图像的数量。我们还考虑了体素间距,即使从同一扫描仪获取,也可能在图像之间变化。这可以通过重新采样到各向同性分辨率来完成:
def resample_img(itk_image, out_spacing=[2.0, 2.0, 2.0], is_label=False):
# Resample images to 2mm spacing with SimpleITK
original_spacing = itk_image.GetSpacing()
original_size = itk_image.GetSize()
out_size = [
int(np.round(original_size[0] * (original_spacing[0] / out_spacing[0]))),
int(np.round(original_size[1] * (original_spacing[1] / out_spacing[1]))),
int(np.round(original_size[2] * (original_spacing[2] / out_spacing[2])))]
resample = sitk.ResampleImageFilter()
resample.SetOutputSpacing(out_spacing)
resample.SetSize(out_size)
resample.SetOutputDirection(itk_image.GetDirection())
resample.SetOutputOrigin(itk_image.GetOrigin())
resample.SetTransform(sitk.Transform())
resample.SetDefaultPixelValue(itk_image.GetPixelIDValue())
if is_label:
resample.SetInterpolator(sitk.sitkNearestNeighbor)
else:
resample.SetInterpolator(sitk.sitkBSpline)
return resample.Execute(itk_image)
# Assume to have some sitk image (itk_image) and label (itk_label)
resampled_sitk_img = resample_img(itk_image, out_spacing=[2.0, 2.0, 2.0], is_label=False)
resampled_sitk_lbl = resample_img(itk_label, out_spacing=[2.0, 2.0, 2.0], is_label=True)
如果需要进一步归一化,我们可以使用医学图像注册包(例如MIRTK等)并将图像注册到相同的空间中,使得图像之间的体素位置彼此对应。分析结构脑MR图像(例如T1加权MR图像)的典型步骤是将训练数据库中的所有图像登记到参考标准,例如平均图集(例如MNI 305)图集)。根据配准方法的自由度,这也可以规范尺寸(仿射配准)或形状(可变形配准)。这两种变体很少使用,因为它们去除了图像中的一些信息(即形状信息或大小信息),这些信息可能对分析很重要(例如,大心脏可能预测心脏病)。
数据增强
通常情况下,可用的数据量是有限的,有些差异没有被涵盖。一些例子包括:
- 软组织器官,有很多正常的形状
- 病理,如癌症病变,在形状和位置上有很大的不同
- 免费的超声波图像,有很多可能的视图。
为了恰当地推广到不可见的测试用例,我们通过模拟我们要对抗的数据的变化来增强训练图像。类似于标准化方法,我们区分强度和空间增强:
强度增强的例子:
- 向训练图像添加噪声可归结为噪声图像
- 添加随机偏移或对比度以处理图像之间的差异
空间增强的例子:
- 在预期对称性的方向上翻转图像张量(例如,在脑部扫描时左/右翻转)
- 随机变形,(例如用于模仿器官形状的差异)
- 沿轴的旋转(例如,用于模拟差异超声视角)
- 随机裁剪
强度和空间增强技术示例
关于扩充和数据I / O的重要说明:根据需要或有用的扩充,某些操作仅在python中可用(例如随机变形),这意味着如果使用使用原始TensorFlow的读取方法(即TFRecords或tf.placeholder),它们需要预先计算并存储到磁盘,从而大大增加了训练数据库的大小。
类平衡
领域专家解释(例如手动分割疾病类别)是在医学图像的监督学习期间的要求。通常,图像级(例如疾病类)或体素级(即分段)标签不能以相同的比率获得,这意味着网络在训练期间将不会从每个类看到相同数量的示例。如果类比率有些相似(例如对于二元分类情况为30/70),则这对准确性没有很大影响。但是,由于大多数损失是整个批次的平均成本,因此网络将首先学会正确预测最常见的类别(例如背景或正常情况,这通常是更多可用的示例)。
训练期间的类不平衡将对罕见现象(例如图像分割中的小病变)产生更大的影响,并且在很大程度上影响测试准确性。
为了避免这种下降,有两种典型的方法来对抗数据集中的类不平衡:
通过采样进行类平衡:在此,我们的目标是在采样期间校正所见示例的频率。这可以通过以下方式来完成:a)从每个类中抽取相等的量,b)过度抽样的过度表示的类或c)过度抽样的频率较低的类。在DLTK中,我们有a)的实现。我们对图像卷中的随机位置进行采样,并考虑一个提取的示例,如果它包含我们正在寻找的类。
通过损失函数进行类平衡:与典型的体素平均损失(例如分类交叉熵,L2等)相比,我们可以a)使用本质上平衡的损失函数(例如平滑骰子损失,这是一种均值所有类别的骰子系数)或b)通过类频率(例如,中频重新加权的交叉熵)对每个预测的损失进行重新加权。
应用亮点示例
通过本博文中提供的所有基本知识,我们现在可以研究使用TensorFlow构建用于医学图像深度学习的完整应用程序。
示例数据集
对于大多数情况,我们使用了IXI脑数据库(http:///ixi-dataset/)。对于图像分割,我们下载了MRBrainS13挑战数据库(http://mrbrains13.isi./),您需要先注册才能下载它。
多通道脑MR图像的图像分割
多序列图像输入,目标标签和预测的Tensorboard可视化
该图像分割应用程序学习在小(N = 5)MRBrainS挑战数据集上从多序列MR图像(T1加权,T1反转恢复和T2 FLAIR)预测脑组织和白质病变。它使用具有残差单位的3D U-Net网络作为特征提取器,并跟踪TensorBoard中每个标签的Dice系数精度。
代码和说明可以在这里找到(https://github.com/DLTK/DLTK/tree/master/examples/applications/MRBrainS13_tissue_segmentation)。
T1加权脑MR图像的年龄回归和性别分类
输入T1加权脑MR图像用于回归和分类示例
采用可扩展3D ResNet架构的两个类似应用程序学习从IXI数据库的T1加权脑MR图像中预测受试者的年龄(回归)或受试者的性别(分类)。这些应用程序之间的主要区别在于损失函数:我们训练回归网络以将年龄预测为具有L2损失的连续变量(预测年龄与实际年龄之间的均方差),我们使用分类交叉 - 熵损失预测性别的类别。
这些应用程序的代码和说明可以在这里找到:分类(https://github.com/DLTK/DLTK/tree/master/examples/applications/IXI_HH_sex_classification_resnet),回归(https://github.com/DLTK/DLTK/tree/master/examples/applications/IXI_HH_age_regression_resnet)。
3T多通道脑MR图像的表示学习
使用深度卷积自动编码器网络测试图像和重建
在这里,我们演示了深度卷积自动编码器架构的使用,这是一种强大的表示学习工具:网络将多序列MR图像作为输入,旨在重建它们。通过这样做,它在其潜在变量中压缩整个训练数据库的信息。训练的权重也可用于转移学习或信息压缩。请注意,重建的图像非常平滑:这可能是由于此应用程序使用L2丢失功能或网络较小以正确编码详细信息。
代码和说明可以在这里找到(https://github.com/DLTK/DLTK/tree/master/examples/applications/IXI_HH_representation_learning_cae)。
T1w脑MR图像上的简单图像超分辨率
图像超分辨率:原始目标图像; 下采样输入图像; 线性上采样图像; 预测超分辨率
单图像超分辨率旨在学习如何从低分辨率输入上采样和重建高分辨率图像。这种简单的实现创建了图像的低分辨率版本,超级分辨率网络学会将图像上采样到其原始分辨率(此处上采样因子为[4,4,4])。此外,我们计算线性上采样版本以显示重建图像的差异。
代码和说明可以在这里找到(https://github.com/DLTK/DLTK/tree/master/examples/applications/IXI_HH_superresolution)。
最后…
我们希望本教程能帮助您轻松深入学习生物医学图像的深度学习。