FSL-VBM流程与说明

更新：2018-03-19
创建：2017-11-06 by 叶辰飞

可选参数已加粗 可选方法见目录

数据准备

原理说明：这里最好从原始结构像生成经过方向矫正以及砍脖子之后的T1结构像。方向矫正会把sagittal plane转成axial plane；砍脖子的效果比较保守，可能会有脖子残留，但不会把脑子砍掉。经测试[^测试]，砍一次脖子和重复砍脖子的效果相差无几，推测不带脖子的原始图像经过此步骤不会被误砍脑子。

对每个subject的原始图像进行方向矫正和砍脖子，代码以nc_daizigang.nii为例：

fslreorient2std nc_daizigang.nii nc_o_daizigang.nii   #o代表reorient
robustfov -i nc_o_daizigang.nii -r nc_co_daizigang.nii    #c代表crop neck

1. 创建路径

原理说明：VBM的数据需求通常要确定两个组别（NC和PA,分别是健康对照和病人组）。每一组包括同类型的T1结构像。每组人数最好多于10人。

新建一个文件夹存放数据

mkdir vbm_foldername
cd vbm_foldername

2. 定义template_list

原理说明：VBM需要创建一个基于用户上传的T1数据的脑模板（specific template），因此首先要选定哪些subject数据将作为脑模板生成的来源。注意为了避免脑模板出现偏差（bias，也就是不能因为某一组人数多，使得脑模板的形态就更偏向某一组），template_list中每组的人数必须相同。这一步在vbm_foldername目录下生成template_list

这一步需要筛选subjects，使得两组的人数一致
因为默认顺序和名字有关，所以名字前缀最好统一为NC_XXX(代表健康对照) 和 PA_XXX(代表病人组)

for g in nc_co_daizigang.nii nc_co_liulixia.nii nc_co_sunchangsheng.nii nc_co_wangshukun.nii nc_co_chenxihe.nii pa_co_songjingxi.nii pa_co_sunfenmiao.nii pa_co_tianyoujin.nii pa_co_wangqingru.nii pa_co_zhaoxiuying.nii  do echo $g >> template_list done

后台测试的时候可以用命令imglob * 来观察subject list是否为NC在前，PA在后
组1默认为NC，组2默认为PA；

QC: 对每一个subject生成平面图,用命令slicesdir `imglob *`保存jpg结果

3. 生成design_matrix

原理说明：design_matrix是一般线性模型（GLM）的setup项,保存了后续统计分析所依赖的数据结构。这一步生成的结果design.mat和design.con在统计分析步骤中会用到。

design_ttest2 design 5 5

这里的5 5分别是两组人的个数，组1组2可以用命令cd ../ \ imglob *查看顺序
注意这里的人数是包含所有人数，不是template_list的人数

4. 剥头皮（skull stripping）

原理说明：剥头皮[1]是很重要的一步，决定了后续皮层灰质分割的效果。因为皮层灰质分割方法假设大脑中只有脑脊液、白质和灰质组成，如果混入了头皮会影响灰质区域的估计。目前没有任何剥头皮算法能在所有数据都表现出高性能，因此这里提供一些可选项。生成的结果是剥头皮后的脑图像XXX_struc_brain.nii。

4.1 方法1：

fslvbm_1_bet -b

这一步是FSL-BET的内置算法，速度最快，但效果一般，默认推荐（FSL内置）。

4.2 方法2：

fslvbm_1_bet -N

也是FSL-BET内置算法，但是不推荐，耗时长。-N适用于输入图像带脖子的情况(前序步骤robustfov的结果一般也会残留少部分脖子，可视第2步QC的结果决定是否采用剥头皮方法2)

4.3 方法3：

mkdir struc
cp -r nc_co_daizigang.nii struc/nc_co_daizigang_struc.nii  
bet struc/nc_co_daizigang_struc.nii struc/nc_co_daizigang_struc_brain.nii -R -S -B

也是FSL-BET内置算法，但是不推荐，耗时非常长。-B去除脖子，-R是多次重复bet，-S去除眼球等软组织；可视第2步QC的结果决定是否采用剥头皮方法3。这里以nc_co_daizigang.nii为例

4.4 方法4：

mkdir struc
cp -r nc_co_daizigang.nii struc/nc_co_daizigang_struc.nii    
~/MindsGo/BrainSuite17a/bin/bse -i struc/nc_co_daizigang_struc.nii -o struc/nc_co_daizigang_struc_brain.nii --auto --trim -p

用BrainSuite的BSE工具剥头皮。这里以nc_co_daizigang.nii为例

4.5 方法5：

mkdir struc
cp -r nc_co_daizigang.nii struc/nc_co_daizigang_struc.nii    
antsBrainExtraction.sh -d 3 -a struc/nc_co_daizigang_struc.nii -e /mnt/hgfs/ShareFolder/Chenfei/Oasis/T_template0.nii.gz -m /mnt/hgfs/ShareFolder/Chenfei/Oasis/T_template0_BrainCerebellumProbabilityMask.nii.gz -o struc/daizigang -f /mnt/hgfs/ShareFolder/Chenfei/Oasis/T_template0_BrainCerebellumRegistrationMask.nii.gz
mv struc/daizigangBrainExtractionBrain.nii.gz struc/nc_co_daizigang_struc_brain.nii.gz

用ANTs的剥头皮工具，效果很好，作为第二推荐。这里以nc_co_daizigang.nii为例。

QC: 原则上，这一步生成的剥头皮结果应该只包含大脑灰质、白质和脑脊液。可以接受残留一点点头皮，或者缺失一点点皮层。如果所有人的剥头皮效果都差，需要换方法；如果只是个别差，考虑把差的subject删除。

5. 生成脑模板

原理说明：首先把每个人剥头皮之后的脑子都进行灰质、白质和脑脊液分割；之后把template_list 中每个人的灰质*_struc_GM都仿射变换配准到FSL内置的标准脑模板上（GM ICBM-152 template）；然后把结果平均，左右镜像，得到一个初级阶段脑模板template_GM_init；之后把每个人的灰质*_struc_GM都配准到template_GM_init（线性或者非线性）；再把得到的结果进行平均，左右镜像，得到左右对称的终极脑模板（specific template，2mmx2mmx2mm分辨率）。生成的结果主要有template_4D_GM和template_GM。

5.1 方法1：

fslvbm_2_template -n

非线性配准，耗时较长，但配准效果好，作为默认推荐

5.2 方法2：

fslvbm_2_template -a

线性配准，耗时较短，但配准效果较差

QC: 可以提供template_4D_GM和template_GM给用户检查。template_4D_GM是每个subject配准到template_GM_init之后的4D结果，理论上每个subject配准后的灰质图像应该处在MNI空间中心位置，允许个体间少量差异存在。template_GM是终极脑模板（specific template），同样也应该处在MNI空间中心位置。

生成的终极脑模板如：

6. 个体配准到模板

原理说明：这一步要把subject list的每个人灰质图配准到上一步生成的脑模板（specific template）上，同时平滑，结果生成在stats文件夹。平滑的目的是增加信噪比。生成的结果主要有GM_mod_merg_s2，GM_mod_merg_s3，GM_mod_merg_s4，GM_mask，fslvbm_tstat1， fslvbm_tstat2。

fslvbm_3_proc

这一步无需任何参数

QC: 可以提供GM_mod_merg_s2,GM_mod_merg_s3和GM_mod_merg_s4给用户检查。末尾的数值代表高斯平滑核的大小（sigma，单位是mm），官方默认推荐s3。原则上，sigma越大，图像越模糊，个体图像的相似性也越高，但是分辨率会更低。因此需要寻找分辨率和信噪比的balance。

7. 统计分析

原理说明：核心思想是置换检验(permutation test)[2]，优点是作为非参数检验所以无需对数据分布做任何假设。理论上，显著性差异更有可能出现在团簇（cluster）中而非孤立的体素（voxel），这里有两种不同的统计方法来利用cluster的信息，均内含了family-wise error矫正（FWE correction）。对于p值阈值，0.95对应corrected p=0.05，0.99对应corrected p=0.01，0.999对应corrected p=0.001。

7.1 方法1：TFCE法[3]（Threshold-Free Cluster Enhancement）

randomise -i GM_mod_merg_s3 -m GM_mask -o fslvbm -d design.mat -t design.con -T -n 5000

方法1作为官方默认推荐。这里以GM_mod_merg_s3为例（用户可选s2/s3/s4），对个体灰质数据置换5000次（用户可自定义）。生成的结果主要有fslvbm_tfce_corrp_tstat1和fslvbm_tfce_corrp_tstat2。这里fslvbm_tfce_corrp_tstat1是最重要的结果，反应了病人组相对于健康对照组的萎缩p值。

以下命令可以在后台测试时快速浏览VBM的P值结果图（如果没有显著结果，可把-b的参数调小，如0.9或0.8）：

fslview $FSLDIR/data/standard/MNI152_T1_2mm fslvbm_tfce_corrp_tstat1 -l Red-Yellow -b 0.95,1

7.2 方法2：cluster-based thresholding法[4]

randomise -i GM_mod_merg_s3 -m GM_mask -o fslvbm -d design.mat -t design.con -c 2.3 -n 5000
fslmaths fslvbm_clustere_corrp_tstat1 -thr 0.95 -bin mask_pcorrected
fslmaths fslvbm_tstat1 -mas mask_pcorrected fslvbm_tstat1_corrected

以下命令可以在后台测试时快速浏览VBM的P值结果图：

fslview $FSLDIR/data/standard/MNI152_T1_2mm fslvbm_tstat1_corrected -l Red-Yellow -b 2.3,4

方法2需要预设基于团簇的矫正阈值（threshold for the cluster-based correction），默认值为2.3（用户可自定义）。这里以GM_mod_merg_s3为例（用户可选s2/s3/s4），首先计算P值图，然后根据P值阈值限制显著性差异结果出现的空间位置。生成结果包括fslvbm_clustere_corrp_tstat1和fslvbm_clustere_corrp_tstat2。这里fslvbm_clustere_corrp_tstat1是最重要的结果，反应了病人组相对于健康对照组的萎缩p值。

8. 文本结果

原理说明：根据上一步得到的显著性区域，通过图谱匹配法，找到其对应的脑区与MNI空间坐标，同时统计出显著cluster的个数。推荐使用哈佛-牛津的皮层图谱和皮层下结构图谱，比较普及。

autoaq -i fslvbm_tfce_corrp_tstat1 -t 0.95 -o report.txt -a "MNI Structural Atlas" 

这里-t的参数给用户提供可选项（0.95, 0.99, 0.999），默认0.95。-a的参数可选，包括FSL中的内置图谱库：

• “Cerebellar Atlas in MNI152 space after normalization with FLIRT”
• “Cerebellar Atlas in MNI152 space after normalization with FNIRT”
• “Harvard-Oxford Cortical Structural Atlas”
• “Harvard-Oxford Subcortical Structural Atlas”
• “JHU ICBM-DTI-81 White-Matter Labels”
• “JHU White-Matter Tractography Atlas”
• “Juelich Histological Atlas”
• “MNI Structural Atlas” （默认)

总结

VBM最后给用户提供的数据下载应至少包括：

文字结果：
report_cortical.txt，report_subcortical.txt (可以在线优化)
参数图结果：
fslvbm_tfce_corrp_tstat1，fslvbm_tfce_corrp_tstat2
fslvbm_clustere_corrp_tstat1，fslvbm_clustere_corrp_tstat2
模板结果：
$FSLDIR/data/standard/MNI152_T1_2mm

参数图应该叠加在模板上显示。如:

前期考虑用截图在线显示；后期考虑加入互动在线显示模块，或者用自己开发的本地软件支持浏览。

附录：QC方法

VBM所有QC相关的基础原始图像文件包括：

• XXX ：每个人的3D T1 （QC标准：方向是否一致，是否用户错传了别的图像，比如错传成T2图）

• XXX_struc_brain ：每个人剥头皮后的3D T1 （QC标准：剥头皮效果是否可接受）

• template_4D_GM：每个subject配准到template_GM_init之后的4D 灰质（QC标准：灰质皮层分割效果是否可接受）

• template_GM：终极脑灰质模板（specific template） （QC标准：生成的灰质模板是否可接受）

• GM_mod_merg_s2,GM_mod_merg_s3和GM_mod_merg_s4： 高斯平滑后的每个subject灰质（QC标准：根据平滑后的效果决定采用哪个图进行后续统计分析。sigma越大，图像分辨率越低，但信噪比越高）

QC平面图生成：

1. 对于3D图像（包括XXX ，XXX_struc_brain，template_GM）：直接使用slicesdir `imglob *`命令生成名为slicesdir的文件夹。可以根据imglob的参数进行文件选择。比如对于XXX_struc_brain的组图QC，命令为slicesdir `imglob *_struc_brain`；对于template_GM这类单个图像QC，命令为slicesdir `imglob template_GM`。最后将每组文件夹的名字更改为对应的3D图像名字。

2. 对于4D图像（包括template_4D_GM ，GM_mod_merg_s2，GM_mod_merg_s3，GM_mod_merg_s4）：对每个4D文件使用命令fslsplit <input> [output_basename],分裂成一串3D文件，每个3D文件的名字默认为[output_basename]0000, [output_basename]0001, [output_basename]0002, ...。顺序通常对应subject list的顺序。按照这个对应关系，把每个3D文件的名字改成subject ID，这样slicesdir之后的html才会显示subject ID。最后将每组文件夹的名字更改为对应的4D图像名字。