最近需要仿真研究T1结构像+RS-fMRI影像融合处理输出目标坐标的路线可行性。就此机会记录下来。
为了完成验证目标处理,首先需要有数据,然后需要准备对应的处理平台和工具箱,进行一系列。那么开始记录~
前言:
为了基于种子点的功能连接分析,找出与种子点的最相关坐标。如下是严老师的培训文档,正是我所需要的。链接如下:
https://d.rnet.co/Course/V3.0EN/4_R-fMRI_Data_Processing_DPARSFA.pdf
https://d.rnet.co/Course/V3.0EN/4_R-fMRI_Data_Processing_DPARSFA.mp4
确实做功能连接是我这次处理的目的,我需要先完成预处理,后才基于种子ROI去做功能连接分析,既可以分析种子ROI和目标脑区的连接情况,也可以分析种子ROI与全脑的功能连接情况。但是可以想象,两个计算量肯定有明显的差异。
万里长征第一步,还是开始尝试吧。
正文:
(1)下载满足要求的开源数据
去OpenNEURO https://openneuro.org/datasets/ds002422/versions/1.1.0
下载数据做仿真验证够够的了,有46个被试,每个被试下面有多个数据记录,完全符合我的需要,有这个做仿真研究真的是太方便了。
(sub-01_T1W.nii就是T1结构像,sub-01_task-rest_bold.nii.gz就是我要的RS-fMRI)
多下载几个数据处理验证没得问题。
也可以去严老师的网站下载数据,不同类别的data可供下载,可以直接下载Demo Data for DPARSF。
https://rfmri.org/content/demonstrational-data-resting-state-fmri
(2)准备数据处理平台
按照如下准备先按照到位,这些是处理和查看处理过程的利器,作为研究核磁影像的人少不了。
安装完Matlab以后,除了MRIcron 是独立软件,其他就都是工具箱。工具箱软件需要先下载到本地,后解压缩放到Matlab 安装路径的toolbox 文件夹下,并打开matlab 界面,从设置路径处导入。这里我是遇到了一个问题,到底是 添加文件夹 还是 添加并包含子文件夹?
那么看下面
澄清一下:如果打算使用DPABI 做预处理,那么这个地方的SPM12 和DPABI的设置路径的方式一定要参考严老师的官网
使用DPABI 一定要去看:http://www.rfmri.org/DPABIErrorHandling
REST 是DPABI 的一个子模块。
(3)实施处理
都说DPABI 是“傻瓜版”的SPM,但是这个用起来还是要花时间去仔细啄米…
先说一下数据摆放路径,非常容易出错,可以先去阅读网站 http://www.rfmri.org/DPABIErrorHandling
如果用DPABI 来做预处理,那数据摆放路径就是讲究,路径关系如果层级发生错误,就会报各种奇怪的错误。
首先新建一个路径文件夹,如202504_MRI是我的文件夹,文件夹下面需要新建两个子文件夹:
如果原始数据都是Dicom 格式,新建T1Raw 和FunRaw 文件夹,并将各被试的T1像拷贝到T1Raw文件夹,RS-fMRI 像拷到FunRaw文件夹。
如果T1是nii格式,RS-fMRI 是.nii.gz 格式就新建T1Img 和FunImg文件夹,并将3个被试的T1像拷贝到T1Img文件夹,RS-fMRI 像拷到FunImg文件夹,我这次的数据就是这一种。
另外,一定要注意,如果使用DPABI进行预处理,这些原始数据的父文件夹是T1Img 和FunImg文件夹,或T1Raw 和FunRaw 文件夹,我曾经将sub当父文件夹出现了难以接受的错误。正确的数据层级关系如下:
当原始数据是.nii+.nii.gz 格式时:
当原始数据是Dicom格式时:
这里插一个处理步骤,因为我看有人讲一个4D.nii.gz 容易报错,需要转化成多个3D ni。
我的fMRI 影像是一个文件4D.nii.gz 文件,然后我通过手动解压以后变成了只变成了一个nii。如果想变成多个3D nii,还是需要工具或者代码实现。下面提供一段可靠的代码4D.nii.gz 转成多个3D nii的matlab:
clc;clear;
fileName = 'sub-03_task-rest_bold.nii.gz';%需解压的fMRI 的原始4D 文件压缩格式文件名,
filepath = 'E:/202504_MRI/';%fMRI 的原始4D 文件路径
gz_file = fullfile(filepath,fileName);tempniiPath = 'E:/202504_MRI/temp/';%初解压缩的nii文件存储路径,临时的,每完成一个数据的解压以后都要建立一个temp
gunzip(gz_file,tempniiPath);%手动初解压缩fileList = dir(fullfile(tempniiPath,'*.nii*'));
outFile = fileList.name;%fMRI 手动初步解压缩后的文件名
inputFilename = fullfile(tempniiPath ,outFile);%解压后的完整文件名outDir = 'E:/202504_MRI/FunImg/sub3/';%解压成3D nii 的存储路径spm('Defaults','fMRI');
spm_jobman('initcfg');
% 创建 batch job
matlabbatch{1}.spm.util.split.vol = {inputFilename};
matlabbatch{1}.spm.util.split.outdir = {outDir}; % 输出目录
matlabbatch{1}.spm.util.split.prefix = 'vol_'; % 输出文件名前缀
spm_jobman('run',matlabbatch);
(2)预处理步骤虽多,但是处理的过程有讲究
在matlab 窗口输入 DPARSFA(是DPARSF advance 版本,支持T1+fMRI 联合分析)后即可进入DPARSFA 交互界面
①模版参数
选中Blank,其实就是空白模版,其他的都有很多默认的参数和处理步骤被选中是有特定的处理,作为新手来说,先从0开始进入吧,选中blank就相当于空白开始,灭有任何的默认选中处理步骤,需要自己去尝试设定探索;
这个地方参数模版,是一些供选择的处理方案,感觉是处理步骤流程模版化,选择不同的模版会在界面上出现对应选中处理步骤,估计也会存在一些差异性的处理流程步骤,如果很清楚自己的处理目标和步骤,就可以直接选中适合的一个,后面只要把数据扫描关联的参数填下就可以run了。
但是这个地方要结合自己的处理目标去分析。
如果只需要做预处理,那只能选择blank,其他的选项中不管是在MNI空间还是original空间 都是默认选中了预处理+后处理。
②工作路径
这个工作路径需要选中建立的文件夹,是一个总的大路径,不要出现中文字符;
③扫描参数
Time Points:采集的时间点数,应该是个整数,预估100~300,具体需要查看数据
TR:是重复时间,这个也是扫描参数,一般这两个参数需要基于fMRI 数据去查看,如果已知就直接填写
这两个参数必须要准确,如果不知道,可以通过软件查看,实在不行就空着,软件会在处理的过程中读取到。
建议在matlab 中通过命令得到点数,这个是.nii.gz 的获取方式
nii=niftiinfo("my-rs-fMRI.nii.gz");
TimePoints = nii.ImageSize(4);
disp(['时间点总数:',num2str(TimePoints)]);
dicom也有其获取参数的方式,比如先转成nii查看,或者用MRIcron/MRIcronGL/3D slicer/FSL /等软件打开查看。
TR(s)这个参数重要,就是重复时间,这个需要基于扫描json文件去查看,也可以用别的软件工具去查看。
④就是将dicom转成nifti格式,由于我的输入就不是dicom,这一步完全不用选。
上面没有标出来,Apply Mats,如果输入是dicom 格式数据,那就选中吧,如果是我这样的nii格式就不选
⑤ 去除因设备开机不稳的时间点,一般会去掉开始的5~10个时间点的数据。这里,我没有去,因为我下载的数据中已经在json中说了手动去掉了3个时间点,我就不再去了,直接填了0,也可以不选中,省的有效数据变短。
⑥时间层矫正,因为功能像是隔层扫描,做这个处理就是要重新把时间对齐,必须要选中;
⑦功能像的扫描层数,就是ImageSize(3),我下载的数据是64X64X36X200,我这里就是36层
后面的Slice Order:扫描次序,我填的[2:2:36,1:2:35](我看了我下载的json文件,我的是隔层扫描,先扫偶数层,后扫奇数层)
⑧参考层,一般用正中间的扫描时刻对应的层作为参考层,也有TR/2最接近的时刻处的层作为参考
我的扫描时间序列
“SliceTiming”: [1.5375,0,1.6225,0.085,1.7075,0.1725,1.7925,0.2575,1.8775,0.3425,1.9625,0.4275,2.05,0.5125,2.135,0.5975,2.22,0.6825,2.305,0.7675, 2.39,0.855,2.475,0.94,2.56,1.025,2.645,1.11,2.7325,1.195,2.8175,1.28,2.9025,1.365,2.9875,1.4525 ],
可以看出0时刻在第2帧,而且偶数层比奇数层的时刻小,说明先偶数后奇数。
⑨ 表明数据开始处理的路径这样可以附加到工作路径文件夹。
⑩ 运行时间层矫正
看看做完时间层矫正后的结果,可以看到matlab 窗口中输出说明我的参数设置正确了,正确完整了时间层矫正,而且还是3个sub都处理了
头动矫正
Realign:必须要选中
Voxel-Specific Head Motion :因为我计算的目标是坐标靶点,对空间定位和连接精度非常敏感,我这里也需要选中,相当于是需要做一个更高级一些些的头动矫正。普通的核磁处理可以不选中。有人说这一步时间长作用不大,我还是先处理看看。
根据上面,我处理的输入是要基于时间层矫正的结果,所以我界面的这个输入文件夹是时间层矫正的输出
处理完以后输出又多了更多文件夹
影像重定向及分割配准
(11 Reorient Fun : 处理fMRI,就是是的头的朝向变正(因为很多情况下存在数据歪斜)
手动调整功能像歪斜,此处需要手动确立AC点。
(12 AutoMask :处理fMRI 自动掩膜生成,就是基于T1结构像,去掉不是大脑区域的部分,生产掩膜图像;
这个步骤其实是对头动矫正以后的fMRI 影像的平均图像进行处理的,可以得到平均图像去除颅骨以后模版mask。
(13 Crop T1: 就是去掉T1中的无关区域,比如脖子,背景、空气、头皮,从而减少图像中的无效空白部分;
(14 Reorient T1: 处理结构像数据,结构像图像重新定向,矫正朝向和歪斜
需要在弹窗中手动标记AC 点。
(15 Bet :T1 结构像处理,剔除颅骨,提取大脑区域
(16 T1 coreg to Fun: 对齐结构像和功能像,用coregistration 将fMRI 与T1 进行对齐,其实就是将T1结构像配准到fMRI 空间,从而实现两者空间对齐。
这个步骤完成了不少内容,首先讲去掉颅骨的的T1,批准到了去掉头皮的fMRI 上,并得到了中间计算参数,可以讲带颅骨的T1配准大片fMRI 空间。
(17 segment :分割,将T1结构像分割成灰质、白质、脑脊液,segment 是传统的分割方法,New segment+DARTEL 是一种新的分割方法(如果T1结构像的质量很好建议选择这个方法)
这个处理步骤,其实是先将已配准到fMRI 平均影像的T1 又配准到标准空间,这样得到了中间形变参数。同时对标准空间下的T1进行分割处理,分割出脑灰质、脑白质、脑脊液等掩膜图像。
(18 分割处理的模版选择,East Asian 是东亚模版比较适合处理东方人的数据(基于东亚数据得到),European是欧美数据模版(基于欧美数据得到,其实就是MNI)
我这个数据下载的,查阅以后还是先选择欧美模版。
注意11和14都是需要手动调整数据朝向和歪斜的,新手还无从下手,暂时先不选。另外如果数据本来就很正,AC-PC都有差不多水平了,而且也没有左右歪斜就要不选。如果要矫正就要研究其选中后的弹窗页面的调整,等我学会了,我再来补充一下。
11和14选中以后都需要手动去差点参考点,会弹窗让手动操作选择。
先看看AC-PC:
上面图示中给出了AC和PC应该落点的位置,可见AC和PC在上下方向几乎是等高的。
在弹出的页面,手动标记参考时,需要讲AC作为参考点,先理解了AC和PC的结构特点就好落点了。
这一步的输出:
下面用软件逐一打开查阅
T1ImgBet 就是T1结构像做了颅骨剔除以后的大脑区域的图像。
如下图可以看出sub1剔除的的比较差,在小脑下还包含了一些其他的组织,sub2有少许杂质组织,sub3剔除的最好,如果演示我肯定用sub3被试的处理结果去演示。
T1ImgCoreg:
这个是T1配准到fMRI 以后的结果,我也打开看看影像的变化
很遗憾,这个文件夹中的数据跟T1Img 中的数据大小和空间方位还有可视化一模一样,说明T1 Coreg to Fun 是不成功了,为啥不成功,我要重新跟踪步骤检查。
这一步要重新到退处理检查了,这一次不能跳过任何的弹窗和提示。
我查阅资料科普一番以后,发现T1 Coreg to Fun,就是T1结构像配准到去头动后的fMRI 的平均图像中,因为我的fMRI 图像是没有做颅骨剥离的,所以此处T 1 Coreg to Fun 输出的应该还是带颅骨的图。跟我原来的理解有出入,先记录一下。
T1ImgNewSegment 文件家中放的是分割的脑灰质、白质、脑脊液文件。
再看看预处理的后半部分,
这部分主要是4个功能:
①用来处理完成头动矫正的fMRI 影像的,就是基于头动还有基线去做回归,去噪
②是标准化
也是用来处理fMRI影像的,在T1配准到fMRI 空间后,其具备了与fMRI 影像同坐标系空间特点,后 在segment 时其又配准到了标准空间,其转换到标准空间的转换参数也代表了fMRI 到标准空间的参数,基于此对fMRI 影像完成标准化。与此同时还做了基于体素大小的重采样;
③空间平滑处理,还是对fMRI 影像处理,FWHM是高斯窗。如果要计算局部一致性REHo,则平滑处理应该放在计算REHO之后。
④带通滤波,两个是一个意思,只是因为滤波既可以先做,也可以后做。因为如果要计算ALFF 和fALFF 就要后做,选后面那个。
因为我只在意最后用来分析计算功能连接,没有计算ALFF/fALFF和分析局部一致性,预处理部分我都是严格按照上面的步骤来的。
感觉计算这两个参数的话,滤波和平滑就要选择后面的了,不能早早选中。
理解了所有的预处理以后确实没怎么报错,只是这个单纯的处理框架看起来很离散,每个处理步骤都有输入和输出还有做的目的都需要分析。
分析三天也才get 到一条初步的处理流程图吧,一直都在边处理边理顺边验证才得到的一张图,也许这条路到最后发现有冗余,但是真的是要吃透才能简化和完善。可以说预处理和准备分析工作占据了80%,后面的功能连接好像没有那么多步骤。以后在具体实施中去应该吧。功能连接分析的路估计就是这么走了。真的挺长的,后面还要好好消化吸收。