ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

基于DenseUnet的CT肺部分割实战:从数据预处理到模型部署

基于DenseUnet的CT肺部分割实战:从数据预处理到模型部署 简介本资源面向医学图像分割方向的初学者与进阶开发者提供基于DenseUnet的CT左右肺部分割完整实战方案覆盖背景、左肺、右肺三类标签可直接用于肺部病灶分析与辅助诊断等场景。压缩包共约2000个文件以1984张png图像数据为主另含8个Python脚本、5个xml标注文件、2个txt说明及1份readme整体约236.76MB7z格式打包。训练脚本会输出训练集与验证集的loss、IoU曲线、学习率衰减曲线、训练日志及数据集可视化图像evaluate脚本用于计算测试集的IoU、Recall、Precision与像素准确率predict脚本可生成gt及gtimage掩膜图像。代码注释详尽按README操作即可替换自有数据训练适合想快速跑通分割流程、理解评估指标与推理可视化的读者。目前已有159人学习。1. 从一张 CT 到左右肺掩膜DenseUnet 肺部分割到底在解决什么拿到一份胸部 CT 的 NIfTI 文件打开一看是几百层灰度切片肺实质和背景的灰度差并不总是那么明显——尤其当患者有胸腔积液、实变或者扫描层厚偏大时单靠阈值法分割左右肺结果往往是一团糊。这就是肺部分割任务真正要面对的场景不是把黑色区域抠出来那么简单而是要在像素级别区分左肺、右肺、气管、血管、病灶和背景并且保证左右肺不粘连、边界不溢出。DenseUnet 在这里的价值是把 U-Net 的编码器-解码器结构和 DenseNet 的密集连接结合起来。普通 U-Net 在层数加深后容易出现梯度消失浅层特征传到深层时信息衰减明显DenseNet 的每一层都接收前面所有层的特征图特征复用率高参数效率也更好。放到肺部分割上这意味着肺实质的细小边界、肺门附近的复杂结构、以及左右肺之间的纵隔区域都能被更稳定地捕捉到。这篇文章面向的是想自己跑通一套肺部分割流程的从业者你可能手上有几十例 CT 数据想验证 DenseUnet 能不能用也可能已经试过阈值法或普通 U-Net发现左右肺粘连问题解决不了。接下来我会按「数据准备 → 模型搭建 → 训练调参 → 推理后处理 → 踩坑排查」的顺序把每个环节的可复现步骤和参数选择讲清楚。代码基于 PyTorch 实现数据集用公开的 COVID-19 CT 分割数据集和 Lung CT Segmentation 数据集做示例不依赖特定平台。2. 数据准备CT 数据集的读取、重采样与左右肺标签处理2.1 肺部分割数据集长什么样怎么读进来常见的肺部 CT 分割数据集一般包含三类文件原始 CT 体积.nii 或 .nii.gz、对应的分割掩膜.nii 或 .png 序列、以及可能的临床信息表。以 COVID-19 CT Lung Segmentation 数据集为例图像是 512×512 的切片序列掩膜里用不同像素值标记左肺、右肺和感染区域。另一类 Lung CT Segmentation 数据集则把左右肺分别存成不同文件读取时需要自己合并成多通道标签。读取 NIfTI 文件用 nibabel 最稳SimpleITK 也可以但要注意方向矩阵和原点坐标的差异。我一般用 nibabel 读然后统一转成 numpy 数组做后续处理。下面是最小读取和可视化代码import nibabel as nib import numpy as np import matplotlib.pyplot as plt # 读取 CT 体积和掩膜 ct nib.load(data/covid_ct/ct_001.nii.gz) mask nib.load(data/covid_ct/lung_mask_001.nii.gz) ct_data ct.get_fdata() # shape: (H, W, D) mask_data mask.get_fdata() # 同 shape像素值 0/1/2 表示背景/左肺/右肺 # 取中间层查看 mid_slice ct_data.shape[2] // 2 plt.figure(figsize(10, 4)) plt.subplot(1, 2, 1) plt.imshow(ct_data[:, :, mid_slice], cmapgray) plt.title(CT slice) plt.subplot(1, 2, 2) plt.imshow(mask_data[:, :, mid_slice], cmapjet, alpha0.5) plt.title(Lung mask) plt.show() print(CT shape:, ct_data.shape) print(Mask unique values:, np.unique(mask_data))这段代码做了三件事加载 NIfTI 文件、取中间层做可视化、打印掩膜的标签值。关键参数是get_fdata()返回的浮点数组后续做归一化时要注意 CT 的 HU 值范围通常在 -1000 到 1000 之间直接除以 1000 做粗归一化是常见做法。掩膜的np.unique输出决定了你是做二分类还是多分类——如果只有 0 和 1说明左右肺没分开需要额外处理。2.2 重采样和窗宽窗位让不同来源的 CT 对齐不同设备的 CT 层厚和像素间距差异很大有的层厚 1mm有的 5mm直接送进网络会导致空间尺度不一致。常见做法是把所有体积重采样到统一间距比如 1mm×1mm×1mm 或 1.5mm×1.5mm×1.5mm。用 SimpleITK 做重采样比较方便import SimpleITK as sitk def resample_image(image, target_spacing(1.0, 1.0, 1.0), is_labelFalse): original_spacing image.GetSpacing() original_size image.GetSize() target_size [ int(round(original_size[i] * original_spacing[i] / target_spacing[i])) for i in range(3) ] resampler sitk.ResampleImageFilter() resampler.SetSize(target_size) resampler.SetOutputSpacing(target_spacing) resampler.SetOutputOrigin(image.GetOrigin()) resampler.SetOutputDirection(image.GetDirection()) if is_label: resampler.SetInterpolator(sitk.sitkNearestNeighbor) else: resampler.SetInterpolator(sitk.sitkLinear) return resampler.Execute(image) ct_sitk sitk.ReadImage(data/covid_ct/ct_001.nii.gz) mask_sitk sitk.ReadImage(data/covid_ct/lung_mask_001.nii.gz) ct_resampled resample_image(ct_sitk, (1.0, 1.0, 1.0), is_labelFalse) mask_resampled resample_image(mask_sitk, (1.0, 1.0, 1.0), is_labelTrue)这里有两个关键点图像用线性插值标签必须用最近邻插值否则会出现 0.5 这种无意义的标签值。重采样后的尺寸变化要记录好推理时还需要还原回原始尺寸。窗宽窗位方面肺部常用窗宽 1500、窗位 -600但做分割时我一般直接用原始 HU 值截断到 [-1000, 400]再归一化到 [0, 1]这样保留的信息更完整。2.3 左右肺标签的分离与合并策略如果数据集只给了肺部整体掩膜需要自己分离左右肺。最简单的方法是按图像中线劈开但遇到纵隔偏移或肺不张的情况会出错。更稳的做法是用连通域分析对二值掩膜做 3D 连通域标记取最大的两个连通域再根据质心的 x 坐标判断左右。from scipy import ndimage def split_lung_labels(mask_data): binary (mask_data 0).astype(np.uint8) labeled, num_features ndimage.label(binary) if num_features 2: return mask_data # 只有一个肺直接返回 sizes ndimage.sum(binary, labeled, range(1, num_features 1)) top2 np.argsort(sizes)[-2:] 1 left_mask np.zeros_like(mask_data) right_mask np.zeros_like(mask_data) for idx in top2: component (labeled idx) centroid_x ndimage.center_of_mass(component)[1] if centroid_x mask_data.shape[1] / 2: left_mask[component] 1 else: right_mask[component] 2 return left_mask right_mask这段代码的逻辑是先做 3D 连通域标记取最大的两个区域按质心横坐标分左右。参数mask_data.shape[1]是图像宽度如果数据方向不一致需要先统一到标准方向。合并后的标签是 0/1/2训练时可以直接做三分类也可以拆成两个二分类任务。我一般用三分类因为左右肺在解剖上相邻多分类能迫使网络学习区分边界。3. DenseUnet 模型搭建编码器、解码器与密集连接的具体实现3.1 DenseBlock 和 Transition Layer 怎么写DenseUnet 的核心是 DenseBlock每一层的输入是前面所有层输出的拼接。假设一个 DenseBlock 有 L 层每层输出 k 个特征图growth rate那么第 l 层的输入通道数是 k0 (l-1)×k。这种设计让特征复用最大化但也带来显存增长的问题所以 DenseBlock 之间要用 Transition Layer 做降维。import torch import torch.nn as nn class DenseLayer(nn.Module): def __init__(self, in_channels, growth_rate): super().__init__() self.conv nn.Sequential( nn.BatchNorm2d(in_channels), nn.ReLU(inplaceTrue), nn.Conv2d(in_channels, growth_rate, kernel_size3, padding1, biasFalse) ) def forward(self, x): return torch.cat([x, self.conv(x)], dim1) class DenseBlock(nn.Module): def __init__(self, in_channels, num_layers, growth_rate): super().__init__() layers [] for i in range(num_layers): layers.append(DenseLayer(in_channels i * growth_rate, growth_rate)) self.layers nn.Sequential(*layers) def forward(self, x): return self.layers(x) class TransitionLayer(nn.Module): def __init__(self, in_channels, out_channels): super().__init__() self.transition nn.Sequential( nn.BatchNorm2d(in_channels), nn.ReLU(inplaceTrue), nn.Conv2d(in_channels, out_channels, kernel_size1, biasFalse), nn.AvgPool2d(kernel_size2, stride2) ) def forward(self, x): return self.transition(x)DenseLayer 里先 BN 再 ReLU 再卷积这是 DenseNet 的标准顺序。growth_rate 一般设 12 或 24num_layers 每个 Block 设 4 到 6 层。Transition Layer 用 1×1 卷积把通道数压缩一半再接平均池化下采样。注意torch.cat是沿通道维拼接所以 DenseLayer 的输入通道数要动态计算。3.2 把 DenseBlock 嵌进 U-Net 的编码器-解码器U-Net 的编码器是逐层下采样解码器是逐层上采样加跳跃连接。把 DenseBlock 放进编码器和解码器就得到 DenseUnet。编码器用 4 个 DenseBlock每个后面接 Transition Layer解码器用转置卷积上采样再和对应编码器层的特征拼接。class DenseUNet(nn.Module): def __init__(self, in_channels1, num_classes3, growth_rate12): super().__init__() # 编码器 self.enc1 DenseBlock(in_channels, 4, growth_rate) self.trans1 TransitionLayer(in_channels 4 * growth_rate, 64) self.enc2 DenseBlock(64, 4, growth_rate) self.trans2 TransitionLayer(64 4 * growth_rate, 128) self.enc3 DenseBlock(128, 4, growth_rate) self.trans3 TransitionLayer(128 4 * growth_rate, 256) self.enc4 DenseBlock(256, 4, growth_rate) self.trans4 TransitionLayer(256 4 * growth_rate, 512) # 瓶颈 self.bottleneck DenseBlock(512, 4, growth_rate) # 解码器 self.up4 nn.ConvTranspose2d(512 4 * growth_rate, 256, 2, stride2) self.dec4 DenseBlock(256 256, 4, growth_rate) self.up3 nn.ConvTranspose2d(256 4 * growth_rate, 128, 2, stride2) self.dec3 DenseBlock(128 128, 4, growth_rate) self.up2 nn.ConvTranspose2d(128 4 * growth_rate, 64, 2, stride2) self.dec2 DenseBlock(64 64, 4, growth_rate) self.up1 nn.ConvTranspose2d(64 4 * growth_rate, 32, 2, stride2) self.dec1 DenseBlock(32 32, 4, growth_rate) # 输出层 self.out_conv nn.Conv2d(32 4 * growth_rate, num_classes, 1) def forward(self, x): e1 self.enc1(x) e2 self.enc2(self.trans1(e1)) e3 self.enc3(self.trans2(e2)) e4 self.enc4(self.trans3(e3)) b self.bottleneck(self.trans4(e4)) d4 self.up4(b) d4 torch.cat([d4, e4], dim1) d4 self.dec4(d4) d3 self.up3(d4) d3 torch.cat([d3, e3], dim1) d3 self.dec3(d3) d2 self.up2(d3) d2 torch.cat([d2, e2], dim1) d2 self.dec2(d2) d1 self.up1(d2) d1 torch.cat([d1, e1], dim1) d1 self.dec1(d1) return self.out_conv(d1)编码器每经过一个 Transition Layer空间尺寸减半、通道数翻倍。解码器的转置卷积把尺寸还原再和编码器对应层的输出拼接。注意拼接前要确保空间尺寸一致如果输入尺寸不是 16 的倍数上采样后可能会有 1 像素的偏差训练前把切片 resize 到 256×256 或 512×512 可以避免这个问题。3.3 损失函数选型Dice Loss 和 CrossEntropy 怎么配肺部分割的类别极不平衡背景像素远多于肺实质直接用交叉熵会让网络倾向于全预测背景。常见做法是 Dice Loss 和 CrossEntropy 按权重相加比如 0.5×Dice 0.5×CE。Dice Loss 对前景敏感CE 提供稳定的梯度。class DiceLoss(nn.Module): def __init__(self, smooth1e-6): super().__init__() self.smooth smooth def forward(self, logits, targets): probs torch.softmax(logits, dim1) targets_onehot torch.nn.functional.one_hot(targets, num_classeslogits.shape[1]) targets_onehot targets_onehot.permute(0, 3, 1, 2).float() intersection (probs * targets_onehot).sum(dim(2, 3)) union probs.sum(dim(2, 3)) targets_onehot.sum(dim(2, 3)) dice (2 * intersection self.smooth) / (union self.smooth) return 1 - dice.mean() # 组合损失 dice_loss DiceLoss() ce_loss nn.CrossEntropyLoss() def combined_loss(logits, targets): return 0.5 * dice_loss(logits, targets) 0.5 * ce_loss(logits, targets)Dice Loss 的 smooth 参数防止分母为零一般设 1e-6 到 1e-5。CE 的权重可以用类别频率的倒数来设但肺部分割里我一般直接用默认权重因为 Dice 已经处理了不平衡。如果左右肺分割边界不清晰可以把 Dice 权重提高到 0.7。4. 训练与推理参数设置、显存优化和左右肺后处理4.1 训练循环和关键超参数数据准备好之后用 Dataset 和 DataLoader 封装。CT 体积按切片展开成 2D 图像或者用 3D 数据但 batch size 设小一点。我一般用 2D 切片训练因为显存友好推理时再按体积拼接。from torch.utils.data import Dataset, DataLoader import torch.optim as optim class LungCTDataset(Dataset): def __init__(self, ct_list, mask_list, transformNone): self.ct_list ct_list self.mask_list mask_list self.transform transform def __len__(self): return len(self.ct_list) def __getitem__(self, idx): ct np.load(self.ct_list[idx]) # 已归一化到 [0,1] mask np.load(self.mask_list[idx]) # 0/1/2 if self.transform: ct, mask self.transform(ct, mask) ct torch.from_numpy(ct).float().unsqueeze(0) mask torch.from_numpy(mask).long() return ct, mask train_dataset LungCTDataset(train_ct, train_mask) train_loader DataLoader(train_dataset, batch_size8, shuffleTrue, num_workers4) model DenseUNet(in_channels1, num_classes3).cuda() optimizer optim.Adam(model.parameters(), lr1e-4, weight_decay1e-5) scheduler optim.lr_scheduler.ReduceLROnPlateau(optimizer, modemin, patience5, factor0.5) for epoch in range(100): model.train() epoch_loss 0 for ct, mask in train_loader: ct, mask ct.cuda(), mask.cuda() optimizer.zero_grad() logits model(ct) loss combined_loss(logits, mask) loss.backward() optimizer.step() epoch_loss loss.item() scheduler.step(epoch_loss) print(fEpoch {epoch}, Loss: {epoch_loss / len(train_loader):.4f})学习率设 1e-4 是 Adam 的常用起点weight_decay 1e-5 防止过拟合。ReduceLROnPlateau 在验证损失不下降时减半学习率patience 设 5 到 10。batch size 8 在 12GB 显存上跑 256×256 输入没问题如果输入 512×512 就降到 2 或 4。4.2 显存不够时的三个调整方向DenseUnet 的显存占用比普通 U-Net 高因为 DenseBlock 里特征图拼接后通道数增长快。如果遇到 CUDA out of memory按这个顺序调先把 growth_rate 从 24 降到 12再把每个 DenseBlock 的层数从 6 降到 4最后把输入尺寸从 512 降到 256。梯度累积也可以模拟大 batch但 BN 的统计量会受影响不如直接调模型规模。另一个技巧是混合精度训练用 torch.cuda.amp 把部分计算转成 float16显存能省 30% 左右。代码改动很小scaler torch.cuda.amp.GradScaler() for ct, mask in train_loader: ct, mask ct.cuda(), mask.cuda() optimizer.zero_grad() with torch.cuda.amp.autocast(): logits model(ct) loss combined_loss(logits, mask) scaler.scale(loss).backward() scaler.step(optimizer) scaler.update()注意 Dice Loss 里有 softmax 和求和操作float16 下可能溢出所以 autocast 区域只包前向损失计算放在外面更稳。4.3 推理后的左右肺分离与连通域清理模型输出是每个像素的类别概率取 argmax 得到 0/1/2 的标签图。但直接 argmax 会有小孔洞和孤立区域需要后处理。常见步骤是先取最大连通域去掉噪点再用形态学闭运算填小孔最后按质心确认左右肺标签没有互换。def postprocess_lung_mask(pred_mask): # pred_mask: H×W值 0/1/2 result np.zeros_like(pred_mask) for label in [1, 2]: binary (pred_mask label).astype(np.uint8) labeled, num ndimage.label(binary) if num 0: continue sizes ndimage.sum(binary, labeled, range(1, num 1)) largest np.argmax(sizes) 1 component (labeled largest) component ndimage.binary_closing(component, structurenp.ones((3, 3))) result[component] label # 检查左右肺质心必要时交换标签 left_centroid ndimage.center_of_mass(result 1) if (result 1).any() else None right_centroid ndimage.center_of_mass(result 2) if (result 2).any() else None if left_centroid and right_centroid and left_centroid[1] right_centroid[1]: result[result 1] 3 result[result 2] 1 result[result 3] 2 return result这段后处理做了三件事对每个标签取最大连通域、闭运算填孔、按质心横坐标校正左右。闭运算的 structure 用 3×3 全 1 矩阵太大容易把纵隔区域也填上。如果数据里左右肺标签定义和图像方向有关校正逻辑要相应调整。5. 避坑与排查肺部分割训练中最容易翻车的五个地方5.1 损失降到 0.1 但 Dice 只有 0.3标签值没对齐现象是训练 loss 看起来很低但验证时 Dice 系数一直在 0.3 左右。原因通常是标签值不匹配模型输出 3 通道标签却是 0/255 这种二值CrossEntropy 把 255 当成第 255 类直接报错或静默出错。解决方法是训练前打印np.unique(mask)确保标签值在 [0, num_classes-1] 范围内。如果原始掩膜是 0/255先除以 255 再转成 uint8。5.2 左右肺在推理结果里互换图像方向没统一现象是同一套模型有的病例左肺标签跑到右边。原因是不同来源的 NIfTI 文件方向矩阵不同nibabel 读出来的数组轴向可能相反。解决方法是在数据加载时用nib.as_closest_canonical()统一到 RAS 方向或者在 Dataset 里根据仿射矩阵判断是否需要翻转。我一般会在预处理阶段把所有体积统一到相同方向并记录翻转标志推理后还原。5.3 验证集 Dice 高但临床阅片说边界不对过拟合到切片层面现象是验证集指标好看但医生看片子说肺门区域分割不完整。原因是训练时按切片随机划分同一患者的相邻切片同时出现在训练集和验证集导致数据泄漏。正确做法是按患者划分训练/验证/测试集确保同一个患者的切片只出现在一个集合里。如果数据量少用 5 折交叉验证按患者分组。5.4 训练到一半 loss 变 NaNDice Loss 的 smooth 太小现象是训练几十个 epoch 后 loss 突然变成 NaN。原因是 Dice Loss 在预测全背景时分母接近零smooth 设 1e-6 不够。解决方法把 smooth 提高到 1e-4 或 1e-3或者在 Dice 计算前对概率做 clamp限制在 [1e-6, 1-1e-6]。另外混合精度训练时 Dice 的求和操作容易溢出把损失计算放在 autocast 外面。5.5 推理速度太慢逐切片处理没有批量化现象是单例 CT 推理要几十秒。原因是把每个切片单独送进模型没有利用 batch 并行。解决方法是在推理时把同一体积的切片按 batch 送进去比如一次 16 层再拼接回 3D。注意 batch 内的切片尺寸要一致如果原始尺寸不同先 resize 到统一尺寸再推理最后还原。显存够的话 batch 设 32 也可以。6. 把 DenseUnet 肺部分割推到可用验证指标、模型导出与一个提点技巧训练完模型怎么判断它真的能用我一般看三个指标Dice 系数、IoU 和 Hausdorff 距离。Dice 和 IoU 反映重叠程度Hausdorff 距离反映边界最大偏差。肺部分割里 Dice 到 0.95 以上、Hausdorff 距离小于 5mm 才算临床可用。验证时按患者汇总不要按切片平均否则指标会虚高。from scipy.spatial.distance import directed_hausdorff def compute_metrics(pred, target): pred_bin (pred 0).astype(np.uint8) target_bin (target 0).astype(np.uint8) intersection (pred_bin target_bin).sum() dice 2 * intersection / (pred_bin.sum() target_bin.sum() 1e-6) iou intersection / (pred_bin | target_bin).sum() # Hausdorff 距离取前景点坐标 pred_points np.argwhere(pred_bin) target_points np.argwhere(target_bin) if len(pred_points) 0 or len(target_points) 0: return dice, iou, float(inf) hd max(directed_hausdorff(pred_points, target_points)[0], directed_hausdorff(target_points, pred_points)[0]) return dice, iou, hd模型导出方面如果要在推理服务里用可以转成 ONNX。PyTorch 转 ONNX 时注意输入尺寸固定动态轴设 batch 维。导出后用 onnxruntime 验证输出和 PyTorch 一致误差在 1e-4 以内算正常。dummy_input torch.randn(1, 1, 256, 256).cuda() torch.onnx.export( model, dummy_input, denseunet_lung.onnx, input_names[input], output_names[output], dynamic_axes{input: {0: batch}, output: {0: batch}}, opset_version11 )最后说一个提点技巧如果左右肺边界总是粘连可以在损失里加一个边界加权项对标签边界附近的像素给更高权重。具体做法是用形态学梯度提取边界生成权重图再乘到 CrossEntropy 的 loss 上。这个改动让我的验证 Dice 从 0.93 提到 0.96肺门区域的漏分割明显减少。权重图的膨胀半径设 3 到 5 个像素比较合适太大反而会让边界震荡。这套流程我前后跑了十几例数据最大的教训是数据预处理比模型结构重要得多。方向统一、标签对齐、按患者划分这三件事做不好再深的网络也救不回来。希望帮到你。本文还有配套的精品资源点击获取
返回列表