已完成脑图像灰白质边界分割，求皮层厚度映射及像素着色方案

问题背景 #

我正在处理脑图像，目标是分离灰质并根据每个点的皮层厚度为其着色（效果类似典型的皮层厚度热力图）。目前已经完成了白质边界与灰质边界的分割，现在需要解决的核心问题是：为每个灰质像素找到最近的白质边界像素，计算欧氏距离，并用这个距离值给灰质区域（包括边界和中间像素）分配厚度颜色。

高效实现思路 #

不要用循环逐个计算像素距离，直接用欧氏距离变换工具，这是处理这类距离计算的标准高效方法，尤其适合3D脑图像场景：

• 先提取白质边界的二值掩码
• 对全图计算到白质边界的欧氏距离图
• 用灰质区域掩码过滤距离图，只保留灰质范围内的厚度值
• 用连续配色可视化厚度分布

完整代码实现 #

import numpy as np
import matplotlib.pyplot as plt
import nibabel as nib
from skimage import filters, morphology
from scipy.ndimage import distance_transform_edt

# 加载原始数据
t1 = nib.load('raw_map1.nii').get_fdata()
t1map = nib.load('thickness_map1.nii').get_fdata()

# 图像预处理
filt_t1 = filters.gaussian(t1, sigma=1)

# 白质边界分割流程
wm = filt_t1 > 75
med_wm = filters.median(wm)
dilw = morphology.binary_dilation(med_wm)
edge_wm = dilw.astype(float) - med_wm
wm_boundary = edge_wm > 0  # 转为二值掩码，白质边界为True

# 灰质区域及边界分割流程
gm = (filt_t1 < 75) & (filt_t1 > 45)
med_gm = filters.median(gm)
dilg = morphology.binary_dilation(med_gm)
edge_gm = dilg.astype(float) - med_gm
dilw2 = morphology.binary_dilation(edge_wm)
fedge_gm = edge_gm.astype(float) - dilw2
fedge_gm2 = fedge_gm > 0
gm_region = med_gm  # 完整灰质区域的二值掩码

# 计算全图到白质边界的欧氏距离图
# 反转掩码：让白质边界为0，其他区域为1，distance_transform_edt会计算每个点到最近0点的距离
distance_map = distance_transform_edt(np.logical_not(wm_boundary))

# 生成灰质专属厚度图：仅保留灰质区域内的距离值，非灰质区域设为0
thickness_map = np.where(gm_region, distance_map, 0)

# 可视化结果
plt.figure(figsize=(12, 6))
plt.subplot(121)
plt.imshow(filt_t1[:,128,:], cmap='gray')
plt.title('预处理后原始脑图')
plt.axis('off')

plt.subplot(122)
im = plt.imshow(thickness_map[:,128,:], cmap='viridis')
plt.title('灰质皮层厚度热力图')
plt.axis('off')
plt.colorbar(im, label='皮层厚度（像素单位）')

plt.tight_layout()
plt.show()

# 可选：保存厚度结果为NIfTI格式
nib.save(nib.Nifti1Image(thickness_map, nib.load('raw_map1.nii').affine), 'cortical_thickness.nii')

关键细节说明 #

• 距离变换优势：scipy.ndimage.distance_transform_edt是优化过的底层实现，比手动循环快百倍以上，完全适配3D脑图像的规模。
• 区域掩码过滤：用np.where把距离图限制在灰质区域内，避免非灰质区域的距离值干扰可视化。
• 配色选择：viridis配色是科学可视化的常用选择，能清晰区分不同厚度的梯度，和目标效果匹配。

内容的提问来源于stack exchange，提问作者flygongaby