MRI脑影像分析从哲学到技术：一文搞懂VBM预处理基本原理（全网最详细解析）

基于体素的形态学方法（voxel-based morphometry, VBM），是分析大脑解剖学（结构）差异最常用方法之一。其通过给大脑volume逐体素打标签（分类）的方式来进行组织分割（segmentation），过程高度自动化，比传统的基于ROI先验假设的分析方式（manual ROI tracing）得到的结果，更加具有 稳定性和可重复性。

VBM分析基于（高分辨率）MRI脑部扫描图像（一般用T1加权图像），涉及的预处理步骤主要包括：空间归一化（spatial normalisation——分割和比较的前提），偏置场校正(bias field correction——降低相同组织的亮度值差异，有利于组织分割)，分割（ segmentation），调制（modulation——把空间归一化过程中产生的变形场（ deformation field）作用到分割结果上，使其中保留原来个体的组织体积）平滑（smoothing ——去噪，弥补分割缺陷，便于统计分析），流程如上图，图源¹。

以下内容，主要源于笔者对NBAlab卢家峰老师的教学视频和PPT的理解性整理和大量的阅读拓展，初学小白一枚，如有错误之处，还望不吝指出。

原创不易，转载请注明出处。

copyright©意疏：https://blog.csdn.net/sinat_35907936/article/details/110087260

空间归一化（spatial normalisation）

为什么要空间归一化——哲学

1、VBM组织分割时，需要使用不同组织的组织概率图（Tissue probability maps，TPM）来作为分割的先验，具体细节见后文。而这些TPMs都是被归一化（配准）到了一个公共模板上了的，所以使用VBM来进行组织分割前提就是，将被试的脑袋归一化到公共模板上。另外，单样本统计分析常常是和一个标准进行比较，而这些标准也是配准到了公共模板的。

2、不同的人，脑袋大小形状都不一样，不同次的扫描，脑袋的摆放位置也不一样，如图1所示，图源¹。而多样本统计分析寻找局部解剖结构变化的显著性，往往只是通过对应位置值相减，再通过假设检验来得到的。

也就是说，如果所有被试的图像相同空间坐标对应的解剖位置不基本一致，即没有一致性，那相减得到的并不是对应组织结构的变化，差异是没有意义的，无法进行进一步组水平统计分析。所以需要让所有被试的脑袋在不失皮质特征差异的前提下，都对齐到公共模板上，以校正大脑整体形状和解剖位置的差异。让不同的被试，不同的扫描之间具有可比性或者一致性。

图1.大脑形状位置差异

归一化到哪——标准模板（template）——技术

公共模板的作用是确定一个公共坐标空间，以让不同的脑袋在该坐标空间下，相同的坐标位置对应的组织结构基本一致。对于特殊的研究对象，比如幼儿，需要通过数据制作专门的公共模板。对于成人，则常常直接使用一些著名的标准模板。

MNI是最常用的标准模板或者标准空间，在NeuroImaging and Surgical Technologies Lab的网站上，我们可以下载最新的一些模板来看看，有人脑也有猴子脑的。

ICBM 152是MNI的标准版，笔者下载了其第六代非线性对称平均脑立体定向配准模型的Nifti格式文件² ，然后用MRIcron查看，结果如图2。由于它通过了152个人脑T1扫描图像的平均处理，所以看起来有点糊，不过我们还是很容易的就可以找到前联合AC和后连合PC点。CAT，FSL，SPM都有用MNI空间作为标准模板。

图2.MNI152标准模板

对于坐标空间，我们知道为了描述三维空间中的某个位置，显然我们需要定义一个三维笛卡尔坐标系（空间），扫描时，Diocm转nii时，配准时都有定义这样的一个坐标系。一个笛卡尔坐标系可以通过坐标原点，轴方向和单位长度唯一确定。MNI标准模板坐标空间Y轴线如下图红线所示，图源³，且其X轴Z轴在正中矢状界面上。AC点在Y轴下方约4mm，并且AC和PC几乎水平对齐，AC点靠近坐标原点。

有的空间坐标原点直接是AC，而且AC-PC line直接是Y轴，如Talairach空间。被试影像在和标准模板配准（归一化）的同时，也将其原始坐标空间转换到了标准模板对应的坐标空间。不同的标准模板对应的坐标空间可能会有差异，但是可以相互转换。

图3.MNI152标准模板坐标空间

copyright©意疏：https://blog.csdn.net/sinat_35907936/article/details/110087260

如何归一化——配准——技术

配准过程中要进行线性和非线性配准。

线性配准就是在做仿射变换，通过乘以一个由贝叶斯后验估计出的变换矩阵来实现，是线性操作。仿射变换包括旋转，平移，缩放和错切¹等，前两者是刚体变换。仿射变换是对图像整体进行操作，所以只能匹配整体的位置和大小，所以配准后的脑袋在形状和组织位置上还有较大差异，如图3左图，图源⁴ 。

非线性配准是在对线性配准后的图像进行扭曲和变形（warping , deforming）操作，将被试图像的脑沟和其他结构尽量对到模板上，使其与模板尽可能的像。操作过程由最优化理论和一些平滑度量控制，目标是最小化平方误差函数¹。非线性配准后不同被试的脑袋形状和组织形状和解剖位置已经变得非常相似——相同空间坐标对应的解剖位置基本一致，如图4右图，图源⁴。

图4.线性配准与非线性配准

最后，如果我们用MNI模板，前面提到MNI标准模板对应的坐标空间原点在AC附近，而AC-PC line靠近Y轴，所以如果我们在采集图像的时候，或者配准前手动校正坐标原点（reorientation）的时候，定义原始坐标空间的原点和轴都应该尽量靠近该标准。那样在进行到这一步时，不仅配准的精度，还有配准的速度都会有所增加。

如果配准结果坏掉了，多半是坐标原点离AC太远或者Y轴严重偏离AC-PC line，此时可以用SPM的Display工具，修改位移量，旋转量来重新定义原始空间的坐标原点。当然，如果对配准精度要求比较高，在配准之前也可以用SPM的Display工具处理一下，以减少在非线性配准的优化过程中落入局部最优解的可能，如图5所示。另外，三个位移量分别对应mm：后面的三个值的相反数，最后填完需要点Reorient保存修改。

图5.SPM的Display工具

Pitch： 表示低头仰头或者点头动作；Yaw： 表示左右转头或者摇头动作；Roll： 表示左右偏头动作；注意到，pitch本意是抛，想想这个动作，就可以记住它表示的方向了。另外，yaw表示的意思是偏航，如果我们走路换方向了，我们的头会不由自主的转动，这样也就记住了它的方位。

另外还有种原因可能导致配准结果坏掉，那就是非线性配准优化过程中过拟合——太追求误差小了，而把某些组织结构给抹掉的，如图5蓝色圈出部分，图源⁵。可以通过加正则项或者惩罚项来减少过拟合，熟悉最优化理论或者机器学习的朋友应该会对这个概念比较熟悉。不过一般的工具包，如SPM，FSL，CAT12等都考虑了这个的，只需点点点就可以了。

图5.非线性配准与过拟合

copyright©意疏：https://blog.csdn.net/sinat_35907936/article/details/110087260

偏置场校正（bias field correction）

什么是偏置场——哲学

扫描设备中磁场不均匀，被试在扫描仪中的位置不同（靠近线圈部分可能更亮），噪声以及许多未知因素，可导致MR图像上同一组织的体素值在不同位置平滑的发生变化，导致这些变化的所有因素统称为偏置场。关于偏置场的定义和解释有很多，这个解释是笔者比较赞同的一个。

为什么要偏置场校正——哲学

一般而言，我们是希望相同的组织体素值是比较均匀的。逐体素的组织分割，是基于每个体素的体素值和该体素位置属于某个组织的先验概率的，先验概率由TPM固定。现在假设在白质区域取两个体素，在TPM上它们属于白质的先验概率都是1，而其中一个的灰度值是0.9，另一个的灰度值却只有0.4，这样前者贝叶斯最大后验分子为0.9，而后者只是0.4，VBM几乎会把两个体素判定为两种组织的，但是它们本属于同一组织，这就会导致严重的分割错误。

如图6中，图源。有些地方白质的亮度已经和灰质很接近了，有些地方的灰质比白质还亮，如果直接对该图进行分割，结果肯定是不理想的。当然实际情况下，不均匀度不会像图中那么大，但处于广泛化的考量，还是会进行偏置场校正。

图6.MR偏置场影响的图像

如何进行偏置场校正——技术

将偏置场估计图像叠加到原图像上，就可以实现偏置场校正，如图7所示，图源⁵。偏置场图像可以通过算法估计得到，如N4BiasFieldCorrection方法。

图7.MR偏置场校正

分割（segmentation）

为什么要进行组织分割——哲学

结构决定功能，不同组织的形态学变化可能与大脑某些功能的变化或者疾病有联系，我们先要获取到不同组织的图像和参数，才能对其进行进一步分析，这需要分割来得到。

如何进行组织分割——技术

最容易理解也是最常用的分割方法之一，是基于高斯混合模型和朴素贝叶斯的。

由中央极限定理，我们可以知道在大样本的条件下，采样值出现的概率是服从高斯分布的。脑影像不同组织可能也包含上万体素，所以它们的值也可以近似认为服从某个高斯分布，如图8所示，图源⁵，当知道每一种组织体素值的均值和方差时，我们可以很容易得到它们的分布，进而很容易用朴素贝叶斯进行分类。

图8.不同脑组织体素值分布

但显然，在脑影像中，我们并不知道某个体素究竟素属于哪个组织，没有办法直接求出其均值和方差，需要参数估计。所有体素都是混在一起的，此时每个体素都服从一个总分布，它不是高斯的，但它可以看做是每个组织服从的高斯分布的叠加，组成所谓的高斯混合模型。其中，EM算法（无监督，聚类）常用于估计混合高斯模型中每一个成分分布的参数，对应脑影像，就是每一个组织的体素值概率分布的参数。

现在假定我们已经估计出每个组织的分布参数了，那就相当于我们已经知道每一种组织出现某个体素值的概率。接下来就可以利用朴素贝叶斯对脑影像中的每一个体素进行打标签工作，判断它们是灰质，白质还是脑脊液。我们知道伟大的贝叶斯公式是下面这样的：

其中，P(B|A)是类条件概率，即每一个组织出现某个体素值的概率，由估计出来的概率分布得到。P(A)是先验概率，即某个空间位置出现某个组织的概率，由TPM得到。P(B)全概率，即某个空间位置出现某体素值的全概率，由全概率公式得到。最后结果就是某个空间位置，出现某个体素值时，它属于某个组织的概率，最大后验贝叶斯将后验概率最大的组织标签，分配给该体素。以上只是基础的方法，在实际的软件中可能会涉及其他很多优化操作。

TPM是标准的先验图谱，它是通过统计得到的，一般来源于众多MNI空间的图像。如图9所示，是灰质，白质，脑脊液和背景的组织概率图，图源⁵。越亮表示该组织在该空间位置出现的可能性越大，即先验概率越大。对于要在原始空间或者自定义的公共空间做分割的情况，可以将该模板逆变换到对应空间。

图8.脑组织概率图

分割后的结果可能是如图9所示这样的。一般只输出灰质和白质的分割结果到两幅图像中，以进行后续的处理。

图9.VBM脑组织分割结果

copyright©意疏：https://blog.csdn.net/sinat_35907936/article/details/110087260

调制（modulation）

为什么要调制——哲学

我们可以知道，空间归一化之后，的确分割与统计分析都变得容易了，因为在相同坐标空间下，不同被试个体相同坐标对应的组织结构基本上是一样的。但是我们会发现一件很诡异的事情，归一化前每一个体素都可以代表一个标准体积单位，但在分割完成后，如果每个体素依旧代表相同的体积单位，那么不同被试个体对应组织的形状体积就变的基本一样了，这样就不用分析体积变化了，反正大家都一样，如果实在有差异，可能也只是配准误差。但不应该是这样的，因为分割后或者说空间归一化后的每一个体素代表的体积已经改变了，或者拓展或者压缩，这由配准过程产生的形变决定，而不再是一个相同的标准体积单位。

我们知道，旋转平移肯定不会使得体积改变，但线性配准中的尺度变换和错切，非线性变换中的扭曲(warp)和变形(deform) 都会使作用后每个体素代表的体积发生变化，而分割后的结果中显然丢失了这种变化的信息。所以我们需要在分割结果中加上这种变化信息，以维持个体组织体积的不变性，这个过程就叫做调制。

如何调制——雅可比行列式——技术

结构形状肯定是不能再改变了，不然空间归一化就白做了。但是体素值可以改变，分割完成后原先的体素值——信号大小，就没有什么作用了，所以笔者觉得分割完成后同一组织的体素值应该是相同的，比如1，它只用于区分区域，没有实际含义。但如果我们将分割后的图像的每一个位置在配准过程中产生的体积变化率乘到其体素值中，对体素值进行缩放，那么由于同一组织在不同位置发生的体积缩放不同，相乘后得到的体素值也将不同，而这些值将含有原先个体的组织体积信息，而它也将有实际意义，体素浓度，代表分割后每个体素含有原空间的体素的个数，或者标准体积单位的个数。

具体化的过程，如图10所示。模板上的某个组织体积和原先的个体组织体积的对比就相当于个体配准后和配准前的组织体积对比。假设原先某个小的组织结构，其体积为1且里面包含体素的体素值也为统一值1，如果模板上该组织结构的体积为2，即配准后该组织体积将变成2，被放大2倍。体积放大两倍的本质是其包含的体素数量变成原来的两倍，为了维持前后体积不变，则新的体素所能代表的标准体积单位缩小为原来的1/2，或者体素浓度变为原来的1/2，被稀释了。将该比率作用到体素值上，则体素值就可以代表体素浓度。

图10.调制具体化

配准过程中体素体积的缩放比例，用雅克比行列式来计算，如下。如果行列式值大于1，则表示此处相较于参考点体积膨胀了，如果行列式值小于1，则表示此处相较于参考点体积压缩了。用它就可以来记录，以归一化后标准空间上每一个体素为参考点时，原先空间相较于标准空间的体积变化率，显然它是和体素浓度变化相同的。

空间的每个位置都会有一个雅克比行列式值，将其与分割后的图像对应位置相乘，就可以得到调制后的分割结果，如图11所示。

图11.VBM分割结果调制

平滑（smoothing）

为什么要平滑——哲学

图像中含有的随机噪声可能带来统计结果的假阳性。

哪怕通过非线性配准配准，图像也不是和模板完全精确对应的。一定程度的平滑，可以使得图像变得模糊，从而弥补配准的不准确。而且世界本来就是马尔科夫的。

用高斯核平滑后的结果，更加服从中央极限定理，有助于统计分析。

如何进行平滑——技术

平滑的本质是将某个体素值与周围的若干体素值加权平均后代替原先的值，使得一个区域内的值更相近，可以用卷积来实现。在用高斯核来进行平滑的时候，模糊程度由其半高全宽（FWHM）决定，一般取8，一个全脑三个维度上就是[8 8 8]。FWHM越大，对应频域就越窄，高频就滤的更严重，图像就会更模糊。二维高斯核如图12所示，图源⁵。

图10.二维高斯核

copyright©意疏：https://blog.csdn.net/sinat_35907936/article/details/110087260

参考

https://www.slideserve.com/thurman/vbm-voxel-based-morphometry ↩︎ ↩︎ ↩︎ ↩︎
http://packages.bic.mni.mcgill.ca/mni-models/icbm152/mni_icbm152_nl_VI_nifti.zip ↩︎
https://mri-q.com/uploads/3/4/5/7/34572113/talairachversusmni.pdf ↩︎
Karl J. Friston et al., Statistical Parametric Mapping, Academic Press, 2006. ↩︎ ↩︎
https://www.slideserve.com/damita/zurich-spm-course-2011-spatial-preprocessing ↩︎ ↩︎ ↩︎ ↩︎ ↩︎