
简介这份资源面向计算机相关专业学生、高校教师及希望入门医学图像分析的开发者提供一套基于Python机器学习的脑PET图像分析与疾病预测完整源码可直接用于课程设计、毕业设计、课程作业或实训实验也适合作为公司项目的二次开发基础。压缩包共6个文件约10KB以4个Python脚本为主分别承担模型训练、疾病预测、图像裁剪与参数配置等核心功能另含1个依赖清单文本和1份项目说明文档便于快速搭建运行环境并理解整体流程。目前已有84人学习关注。读者可借助训练与预测脚本掌握从数据预处理、特征提取到模型调优的完整链路通过裁剪脚本定位感兴趣区域并依据配置文件和依赖清单复现实验遇到问题还可与作者沟通获取远程指导适合作为医学图像分析与机器学习交叉方向的入门实践参考。1. 脑PET图像分析做疾病预测从一套Python源码和数据集能跑出什么脑PET图像分析和疾病预测这个方向真正动手做过的人都知道难点从来不是模型结构本身而是数据怎么读、标签怎么对齐、预处理怎么统一。一套标着「源码数据集」的Python机器学习项目价值不在于它用了多深的网络而在于它把从DICOM或NIfTI读片、切片筛选、特征归一化到分类器训练的整条链路串起来了。脑PET反映的是脑区代谢活性阿尔茨海默、轻度认知障碍这类疾病在图像上表现为特定脑区代谢降低机器学习要做的就是把这个空间模式学出来。这套东西适合两类人一类是刚入门机器学习、想找一个有真实医学图像背景的完整项目练手的人另一类是做医学图像分析、需要一套可复现baseline来对比自己方法的人。下面按「数据怎么进、模型怎么搭、坑在哪、怎么验证」的顺序拆开讲。2. 脑PET数据集的结构与预处理读片、配准、切片筛选2.1 先搞清楚数据集里到底有什么拿到一套脑PET图像数据集第一件事不是写模型是把目录结构和文件格式摸清楚。常见的组织方式有两种一种按受试者分文件夹每个文件夹里是该受试者的一次或多次扫描另一种是已经切好的二维切片按类别分目录。PET图像常见格式是DICOM原始临床格式和NIfTI.nii/.nii.gz分析常用格式。如果数据集里是DICOM你需要用pydicom逐个读取如果是NIfTI用nibabel更直接。先跑一段探查代码把文件数量、维度、体素间距、标签分布打印出来。这一步不做后面全是盲猜。import os import nibabel as nib import numpy as np from collections import Counter data_root ./data # 数据集根目录按实际路径改 # 统计文件格式和标签分布 ext_counter Counter() label_counter Counter() shapes [] for root, dirs, files in os.walk(data_root): for f in files: ext f.split(.)[-1].lower() ext_counter[ext] 1 if ext in (nii, gz): path os.path.join(root, f) try: img nib.load(path) shapes.append(img.shape) # 假设标签在父目录名里如 AD / NC / MCI label os.path.basename(root) label_counter[label] 1 except Exception as e: print(读取失败:, path, e) print(文件格式分布:, ext_counter) print(标签分布:, label_counter) print(图像维度样例:, shapes[:5])这段代码做三件事遍历目录统计扩展名确认数据格式用nibabel尝试加载每个NIfTI文件记录维度从父目录名推断标签并统计类别分布。参数上data_root要改成你自己的路径如果标签存在CSV里而不是目录名里label那行要改成查表。跑完你就能看到类别是否均衡——如果AD和NC数量差三倍以上后面训练必须做重采样或加权。2.2 配准与强度归一化PET预处理的三个必调参数PET图像预处理里配准和归一化是绕不开的。不同受试者的脑在空间位置上不一致直接送进模型等于让模型学位置噪声。常见做法是把所有图像配准到MNI标准空间工具用ANTs或FSL的flirt。如果你不想装这些重型工具Python里可以用nilearn的resample_img做重采样到统一模板。强度归一化更关键。PET的SUV值标准化摄取值有物理意义但很多数据集给的是原始计数不同扫描仪、不同注射剂量下数值范围差异很大。常见做法是除以全脑均值或者做z-score。我一般会先试除以全脑均值因为PET的代谢模式是相对值更有意义。import nibabel as nib import numpy as np from nilearn.image import resample_img def preprocess_pet(path, target_affineNone, target_shapeNone): img nib.load(path) data img.get_fdata().astype(np.float32) # 步骤1重采样到统一空间如果提供了目标affine if target_affine is not None: img resample_img(img, target_affinetarget_affine, target_shapetarget_shape, interpolationcontinuous) data img.get_fdata().astype(np.float32) # 步骤2强度归一化——除以全脑均值 brain_mask data np.percentile(data, 20) # 粗略去掉背景 if brain_mask.sum() 0: mean_val data[brain_mask].mean() if mean_val 0: data data / mean_val # 步骤3裁剪到固定尺寸方便批量训练 data np.clip(data, 0, 5) # 去掉极端离群值 return data # 示例加载并预处理 vol preprocess_pet(./data/subject01/pet.nii.gz) print(预处理后形状:, vol.shape, 数值范围:, vol.min(), vol.max())三个参数要重点调np.percentile(data, 20)这个阈值决定脑掩膜范围PET背景噪声大阈值太低会把噪声算进均值太高会切掉低代谢区np.clip(data, 0, 5)的上限5是经验值归一化后脑区代谢一般在0到3之间超过5的基本是噪声或伪影重采样的target_shape要和你的模型输入对齐常见是(128,128,128)或(96,96,96)太大显存吃不消太小丢细节。2.3 切片筛选不是每一层都值得送进模型三维PET体积直接送进3D CNN对显存要求高很多项目会转成2D切片做分类。但脑PET的切片不是每层都有信息——顶部和底部的切片大部分是背景或颅骨外区域。常见做法是只取中间一定比例的切片或者按脑区代谢总量筛选。我一般会计算每层切片的非零体素占比保留占比超过某个阈值的层。这样能自动去掉空层减少无效计算。def select_slices(vol, min_ratio0.15): 按非零体素占比筛选有效切片 selected [] for i in range(vol.shape[2]): sl vol[:, :, i] nonzero_ratio (sl 0.1).sum() / sl.size if nonzero_ratio min_ratio: selected.append(i) return selected slices select_slices(vol) print(f原始层数: {vol.shape[2]}, 筛选后: {len(slices)})min_ratio0.15这个阈值需要根据你的数据调。如果数据集已经做过脑提取0.1到0.2都合理如果没做脑提取背景占比高可能要降到0.05。筛选后如果层数还很多可以再做等间隔采样比如每隔2层取1层。3. 用Python搭脑PET疾病预测模型从特征工程到分类器3.1 选2D CNN还是3D CNN显存和精度的权衡这是做脑PET分类第一个要做的选型决定。3D CNN能直接学空间三维模式理论上更适合PET这种 volumetric 数据但对显存要求高batch size往往只能开到2到4训练不稳定。2D CNN把每个切片当独立样本可以用大batch训练快但丢失了层间关系。我的经验是如果数据集样本量在几百以内优先用2D CNN加切片聚合比如每个受试者取多张切片预测后投票或平均如果样本量上千且有条件上大显存再考虑3D。常见做法是用预训练的2D backbone如ResNet18做迁移学习因为医学图像数据量通常不够从头训。下面给一个2D CNN的完整训练框架用PyTorch。结构不复杂重点在数据加载和标签对齐。import torch import torch.nn as nn from torch.utils.data import Dataset, DataLoader import numpy as np class PETSliceDataset(Dataset): def __init__(self, volumes, labels, slice_indices): volumes: list of 3D numpy arrays (已预处理) labels: list of int slice_indices: list of list, 每个受试者选中的切片索引 self.samples [] for vol, label, idxs in zip(volumes, labels, slice_indices): for i in idxs: self.samples.append((vol[:, :, i], label)) def __len__(self): return len(self.samples) def __getitem__(self, idx): sl, label self.samples[idx] sl np.expand_dims(sl, 0) # 加通道维 return torch.FloatTensor(sl), torch.LongTensor([label])[0] class SimplePETNet(nn.Module): def __init__(self, num_classes3): super().__init__() self.features nn.Sequential( nn.Conv2d(1, 32, 3, padding1), nn.BatchNorm2d(32), nn.ReLU(), nn.MaxPool2d(2), nn.Conv2d(32, 64, 3, padding1), nn.BatchNorm2d(64), nn.ReLU(), nn.MaxPool2d(2), nn.Conv2d(64, 128, 3, padding1), nn.BatchNorm2d(128), nn.ReLU(), nn.AdaptiveAvgPool2d(1) ) self.classifier nn.Linear(128, num_classes) def forward(self, x): x self.features(x) x x.view(x.size(0), -1) return self.classifier(x) # 训练循环 def train_model(model, loader, epochs20, lr1e-3): optimizer torch.optim.Adam(model.parameters(), lrlr) criterion nn.CrossEntropyLoss() model.train() for epoch in range(epochs): total_loss 0 for x, y in loader: optimizer.zero_grad() out model(x) loss criterion(out, y) loss.backward() optimizer.step() total_loss loss.item() print(fEpoch {epoch1}, Loss: {total_loss/len(loader):.4f})关键参数说明num_classes3对应AD/NC/MCI三类按你的标签数改lr1e-3是Adam的常用起点如果loss震荡就降到1e-4epochs20只是示例实际要看验证集早停。数据加载部分slice_indices来自上一章的筛选结果每个受试者贡献多张切片标签继承受试者标签。3.2 类别不均衡与数据泄漏两个最容易翻车的地方医学数据集类别不均衡是常态。AD和NC可能各几百例MCI只有几十例。直接训练模型会偏向多数类。常见处理有三种加权损失函数、过采样少数类、或者分层采样。我一般先用加权损失简单且不改数据分布。from sklearn.utils.class_weight import compute_class_weight # 假设 labels 是全部受试者的标签列表 class_weights compute_class_weight(balanced, classesnp.unique(labels), ylabels) weights torch.FloatTensor(class_weights) criterion nn.CrossEntropyLoss(weightweights)compute_class_weight的balanced模式会自动按类别频率反比给权重少数类权重高。把weights传给CrossEntropyLoss就行。数据泄漏是另一个血泪坑。同一个受试者可能有多张切片如果随机划分训练集和测试集同一受试者的切片可能同时出现在两边测试精度虚高。正确做法是按受试者划分确保同一个人的所有切片只出现在一个集合里。from sklearn.model_selection import GroupShuffleSplit # groups 是每个切片对应的受试者ID gss GroupShuffleSplit(n_splits1, test_size0.2, random_state42) train_idx, test_idx next(gss.split(samples, labels, groupssubject_ids))groupssubject_ids保证同一受试者的切片不会被分到两边。这个细节不做论文里的精度全是假的。3.3 特征工程路线HOG、灰度共生矩阵与PCA降维如果你不想上深度学习或者数据量太小传统特征工程加分类器也是一条路。脑PET的代谢模式可以用纹理特征描述常用的是灰度共生矩阵GLCM和HOG。提取每个切片的纹理特征拼接后做PCA降维再送SVM或随机森林。from skimage.feature import graycomatrix, graycoprops from sklearn.decomposition import PCA from sklearn.svm import SVC from sklearn.pipeline import Pipeline def extract_glcm_features(slice_2d): # 量化到16级减少计算量 img ((slice_2d - slice_2d.min()) / (slice_2d.max() - slice_2d.min() 1e-8) * 15).astype(np.uint8) glcm graycomatrix(img, distances[1, 3], angles[0, np.pi/4, np.pi/2], levels16, symmetricTrue, normedTrue) feats [] for prop in [contrast, dissimilarity, homogeneity, energy]: feats.extend(graycoprops(glcm, prop).flatten()) return np.array(feats) # 构建pipeline pipeline Pipeline([ (pca, PCA(n_components50)), (svm, SVC(kernelrbf, C1.0, gammascale)) ])distances[1,3]控制共生矩阵的步长步长太小只捕捉局部噪声太大跨脑区levels16是量化级数PET动态范围大量化太细GLCM稀疏太粗丢信息PCA的n_components50是经验值一般保留到解释方差90%以上。SVM的C和gamma用网格搜索调别手拍。4. 脑PET疾病预测的避坑与排查五个真实踩过的坑4.1 坑一图像方向搞反左右脑颠倒现象模型在训练集上精度很高换一批数据就崩可视化发现激活区域和已知病理脑区对不上。原因NIfTI文件的affine矩阵决定了图像方向不同工具保存的方向可能不同。nibabel加载后如果直接取data可能得到左右颠倒或上下颠倒的数组。解决加载后统一用nib.as_closest_canonical(img)把方向转到标准RAS再取data。这一步加在预处理最前面别省。4.2 坑二归一化用了全体数据统计量现象验证集精度远低于训练集且验证集输入数值范围和训练集不一致。原因做z-score时用了整个数据集的均值和标准差测试集信息泄漏到训练过程。解决归一化统计量只能从训练集算然后应用到验证集和测试集。如果做交叉验证每一折重新算。4.3 坑三切片筛选阈值设太高MCI病例被筛没现象MCI类样本在筛选后剩得极少模型完全学不到MCI模式。原因MCI的代谢降低不如AD明显切片非零占比可能偏低统一阈值把MCI的有效层筛掉了。解决按类别分别统计切片占比分布阈值按类别调整或者对MCI类降低阈值。更稳妥的做法是保留所有层让模型自己学权重。4.4 坑四batch size太小导致BN失效现象3D CNN训练时loss剧烈震荡BN层统计量不稳定。原因3D模型显存占用大batch size只能开到2BatchNorm在极小batch下方差估计不准。解决换GroupNorm或InstanceNorm或者用梯度累积模拟大batch。2D路线的话batch size至少开到16。4.5 坑五测试集受试者与训练集重叠现象测试精度95%以上实际部署惨不忍睹。原因按切片随机划分同一受试者的不同切片分到了训练和测试两边。解决按受试者ID做GroupShuffleSplit确保受试者级别隔离。这个坑最隐蔽因为代码不报错指标还很好看。5. 验证脑PET疾病预测模型是否真的学到东西Grad-CAM与置换检验模型精度高不代表学到的是病理特征。脑PET分类模型可能学到的是扫描仪差异、头动伪影、甚至数据集来源的批次效应。要验证模型是否真的关注脑区代谢模式我一般做两件事Grad-CAM可视化和置换检验。Grad-CAM能告诉你模型做决策时看的是图像哪个区域。如果热力图集中在脑室、白质或颅骨外说明模型没学到东西。下面是一个2D CNN的Grad-CAM实现。import torch import torch.nn.functional as F import numpy as np import cv2 def grad_cam(model, input_tensor, target_layer): input_tensor: (1, 1, H, W) features [] grads [] def forward_hook(module, inp, out): features.append(out) def backward_hook(module, grad_in, grad_out): grads.append(grad_out[0]) h1 target_layer.register_forward_hook(forward_hook) h2 target_layer.register_full_backward_hook(backward_hook) model.eval() output model(input_tensor) pred_class output.argmax(dim1).item() model.zero_grad() output[0, pred_class].backward() feat features[0].detach().numpy()[0] # (C, H, W) grad grads[0].detach().numpy()[0] # (C, H, W) weights grad.mean(axis(1, 2)) # (C,) cam np.zeros(feat.shape[1:], dtypenp.float32) for i, w in enumerate(weights): cam w * feat[i] cam np.maximum(cam, 0) cam cv2.resize(cam, (input_tensor.shape[3], input_tensor.shape[2])) cam (cam - cam.min()) / (cam.max() - cam.min() 1e-8) h1.remove() h2.remove() return cam, pred_class # 使用取模型最后一个卷积层 target_layer model.features[6] # 对应第三个Conv2d后的ReLU cam, pred grad_cam(model, sample_input, target_layer) print(预测类别:, pred, 热力图范围:, cam.min(), cam.max())target_layer选最后一个卷积层太浅的热力图太粗糙太深的尺寸太小。cv2.resize把热力图还原到输入尺寸。拿到cam后叠加到原始切片上看激活区域是否落在颞叶、顶叶这些AD典型受累区。置换检验更直接把测试集标签随机打乱重新算精度。如果打乱后精度还在50%以上三分类基线33%说明模型可能靠数据泄漏或批次效应在猜。正常情况打乱后精度应该掉到随机水平附近。import numpy as np from sklearn.metrics import accuracy_score def permutation_test(model, test_loader, n_permutations100): model.eval() all_preds, all_labels [], [] with torch.no_grad(): for x, y in test_loader: out model(x) all_preds.extend(out.argmax(dim1).numpy()) all_labels.extend(y.numpy()) all_preds np.array(all_preds) all_labels np.array(all_labels) real_acc accuracy_score(all_labels, all_preds) perm_accs [] for _ in range(n_permutations): shuffled np.random.permutation(all_labels) perm_accs.append(accuracy_score(shuffled, all_preds)) p_value np.mean(np.array(perm_accs) real_acc) return real_acc, p_value real_acc, p_val permutation_test(model, test_loader) print(f真实精度: {real_acc:.4f}, 置换检验p值: {p_val:.4f})n_permutations100够用p值小于0.05说明模型精度不是随机猜出来的。这个检验花不了多少时间但能帮你判断结果可不可信。最后说个习惯我做完任何医学图像分类项目都会先把Grad-CAM叠加图打印出来看一遍再决定要不要继续调参。如果热力图不对调参调到天荒地老也是白搭。希望帮到你。本文还有配套的精品资源点击获取