心理生理交互作用(PPI)分析批处理(含ROI制作)

引言

先前简单介绍了PPI分析的基本原理,提到了ROI选取的问题,本文就ROI选取以及mask制作的流程进行详细阐述。为便于以后查阅,也将之前的PPI分析批处理code放在一起。

一、ROI选取

在SPM12的示例数据中,针对单个被试的PPI处理,选取ROI“V2”时,是采用的脑激活图与V2的peak点坐标为核心的半径为6mm的球的交集。

%示例数据VOI提取code
matlabbatch{5}.spm.util.voi.name = 'V2';
matlabbatch{5}.spm.util.voi.roi{1}.spm.spmmat = {''};
matlabbatch{5}.spm.util.voi.roi{1}.spm.contrast = 3;
matlabbatch{5}.spm.util.voi.roi{1}.spm.threshdesc = 'FWE';
matlabbatch{5}.spm.util.voi.roi{1}.spm.thresh = 0.05;
matlabbatch{5}.spm.util.voi.roi{1}.spm.extent = 0;
matlabbatch{5}.spm.util.voi.roi{2}.sphere.centre = [15 -78 -9];
matlabbatch{5}.spm.util.voi.roi{2}.sphere.radius = 6;
matlabbatch{5}.spm.util.voi.roi{2}.sphere.move.local.spm = 1;
matlabbatch{5}.spm.util.voi.expression = 'i1 & i2';

不过由于示例数据只是做了单个被试,而我们分析的是两组,具有对照组。所以参考吴士豪的博士论文,我们决定将两组组分析结果的激活图进行重叠,重叠后的区域再和选定的brodmann功能区或者说AAl解剖区域进行重叠,最后得到的mask作为感兴趣区域。

二、ROImask的制作

1.将RW组与SD组正激活、负激活的cluster在DPABI中分别保存为mask,得到文件如下:


1

2.这里演示两组正激活图像重叠code过程:

%-------------------------------
%时间:20220922
%程序员:wishing
%目的:求两张激活图cluster.nii的交集,
%--------------------------------------------
clear,clc;
cd 'D:\desk\speed\data_1\Dp_spm\PPI\ROImask\mask_make';
%重叠RW与SD两组的激活图
%读取两个nii文件信息
W1=spm_vol('rwred_mask.nii');
W2=spm_vol('sdred_mask.nii');
%读取两个nii文件的体积
Y1=spm_read_vols(W1);
Y2=spm_read_vols(W2);
%展示nii文件信息
%imshow(Y1(:,:,30));
%取交集
Y2=fillmissing(Y2,'constant',0);%将Y2中的nan值替换为0,若不存在则不替换
YY=Y1&Y2;
%命名,保存为mask.nii
V=W1;
V.fname='RSred_mask.nii';
spm_write_vol(V,YY);

3.根据经验选取BA9区域,提取MRIcron中Brodmann模板的BA9区
1).打开DPABI的view界面

2).overlay选择MRIcron中的BA模板

1

2

3).设置阈值,设置9就是BA9区,然后保存cluster就能得到BA9的mask

3

4.改变BA9mask文件的维度,使其和SPM处理结果的文件维度相同、
1)这里需要准备一个SPM处理结果的文件,这里我选取了GLM分析后的beta_0009.nii.
2)使用函数处理。

P1='BA9.nii';P2='beta_0009.nii';  
spm_mask(P1,P2);
%这样P1的维度就和P2一样了,而且生成了一个mP2文件,
%注意此函数格式为: FORMAT spm_mask(P1, P2, thresh)
%  P1     - matrix of input image filenames from which
 %        to compute the mask.
%P2     - matrix of input image filenames on which
%           to apply the mask.
 %thresh - optional threshold(s) for defining the mask.
 %The masked images are prepended with the prefix `m'.
%有时需要设置thresh为1

5.再将重叠后的图像与BA9(根据经验选取的模板)重叠,这里与第2步的code类似,注意BA9图像中可能有空值,这时就需要用到Y2=fillmissing(Y2,'constant',0);这个函数。
6.最后得到生成的mask文件。

三、批量生成批处理需要的文件夹

%时间:20220919
%程序员:wishing
%优化后的批量制作文件夹
%--------------------------------------------------------------------------
% 设置循环
clear,clc;
%————————————————————————————————
for i=1:37
p(i) = string(num2str(i, "%02d"));%生成“01”“02”一直到“37”的字符串
end
%--------------------------------------------------------------------------
% 设置文件夹
%————————————————————————————————
for n=1:length(p)
   filename=string(append('D:\desk\speed\data_1\Dp_spm\PPI\RSBA18\RW\sub_1',p(n)));
   mkdir(filename);
end

四、PPI批处理流程code

%-----------------------------------------------------------------------
% Job saved on 19-Sep-2022 11:16:50 by cfg_util (rev $Rev: 7345 $)
% spm SPM - SPM12 (7771)
% cfg_basicio BasicIO - Unknown
%时间:20220919
%程序员:wishing
%目的:实现任务态脑区激活分析后的PPI全流程循环处理,包含VOI提取,PPI变量制作,PPIGLM模型设计与contrast设计

%-----------------------------------------------------------------------
% Initialise SPM
%--------------------------------------------------------------------------
clear,clc;
spm('Defaults','fMRI');
spm_jobman('initcfg');
% spm_get_defaults('cmdline',true);
%--------------------------------------------------------------------------
% 设置循环
%—————————————————————————————————————
for i=1:37
p(i) = string(num2str(i, "%02d"));%生成“01”“02”一直到“37”的字符串
end
for n=1:length(p)
% 设置文件路径
%--------------------------------------------------------------------------
filename_out=char(append('D:\desk\speed\data_1\Dp_spm\PPI\RSBA18\RW\sub_1',p(n)));
filename_spm=char(append('D:\desk\speed\data_1\Dp_spm\RW\FunImg\sub_1',p(n),'\SPM.mat'));
filename_in=char(append('D:\desk\speed\data_1\Dp_spm\RW\FunImg\sub_1',p(n),'\'));
filename_rp=char(append('D:\desk\speed\data_1\Dp_spm\RW\FunImg\sub_1',p(n),'\rp_1',p(n),'_speed_00001.txt'));
file_ppi=char(append('D:\desk\speed\data_1\Dp_spm\RW\FunImg\sub_1',p(n),'\PPI_s-b_RSBA18.mat'));
file_mask='D:\desk\speed\data_1\Dp_spm\PPI\ROImask\mask_make\mask\RSred_ba18_mask.nii,1';
%--------------------------------------------------------------------------
% VOI提取
%--------------------------------------------------------------------------
matlabbatch{1}.spm.util.voi.spmmat = {filename_spm};
matlabbatch{1}.spm.util.voi.adjust = 0;
matlabbatch{1}.spm.util.voi.session = 1;
matlabbatch{1}.spm.util.voi.name = 'RSBA18';
% matlabbatch{1}.spm.util.voi.roi{1}.spm.spmmat = {''};
% matlabbatch{1}.spm.util.voi.roi{1}.spm.contrast = 1;
% matlabbatch{1}.spm.util.voi.roi{1}.spm.threshdesc = 'FWE';
% matlabbatch{1}.spm.util.voi.roi{1}.spm.thresh = 0.05;
% matlabbatch{1}.spm.util.voi.roi{1}.spm.extent = 0;
matlabbatch{1}.spm.util.voi.roi{1}.mask.image = {file_mask};
matlabbatch{1}.spm.util.voi.roi{1}.mask.threshold = 0.5;
matlabbatch{1}.spm.util.voi.expression = 'i1';
%--------------------------------------------------------------------------
% PPI交互作用变量生成
%--------------------------------------------------------------------------
matlabbatch{2}.spm.stats.ppi.spmmat = {filename_spm};
matlabbatch{2}.spm.stats.ppi.type.ppi.voi(1) = cfg_dep('Volume of Interest:  VOI mat File', substruct('.','val', '{}',{1}, '.','val', '{}',{1}, '.','val', '{}',{1}), substruct('.','voimat'));
matlabbatch{2}.spm.stats.ppi.type.ppi.u = [1 1 1
                                           2 1 -1];%这里还不确定,第一列是condition数,第二列是con_name,第三列是权重
matlabbatch{2}.spm.stats.ppi.name = 's-b_RSBA18';
matlabbatch{2}.spm.stats.ppi.disp = 0;
%--------------------------------------------------------------------------
% GLM PPImodel设计
%--------------------------------------------------------------------------
matlabbatch{3}.spm.stats.fmri_spec.dir = {filename_out};
matlabbatch{3}.spm.stats.fmri_spec.timing.units = 'secs';
matlabbatch{3}.spm.stats.fmri_spec.timing.RT = 2;
matlabbatch{3}.spm.stats.fmri_spec.timing.fmri_t = 35;
matlabbatch{3}.spm.stats.fmri_spec.timing.fmri_t0 = 18;
%session
f = spm_select('FPList',filename_in, '^swr.*\.nii$');
matlabbatch{3}.spm.stats.fmri_spec.sess.scans = cellstr(f);
%regressors PPI交互项和头动参数
matlabbatch{3}.spm.stats.fmri_spec.sess.cond = struct('name', {}, 'onset', {}, 'duration', {}, 'tmod', {}, 'pmod', {}, 'orth', {});
matlabbatch{3}.spm.stats.fmri_spec.sess.multi = {''};
matlabbatch{3}.spm.stats.fmri_spec.sess.regress = struct('name', {}, 'val', {});
matlabbatch{3}.spm.stats.fmri_spec.sess.multi_reg = {...
    file_ppi;...
    filename_rp};
matlabbatch{3}.spm.stats.fmri_spec.sess.hpf = 192;
matlabbatch{3}.spm.stats.fmri_spec.fact = struct('name', {}, 'levels', {});
matlabbatch{3}.spm.stats.fmri_spec.bases.hrf.derivs = [0 0];
matlabbatch{3}.spm.stats.fmri_spec.volt = 1;
matlabbatch{3}.spm.stats.fmri_spec.global = 'None';
matlabbatch{3}.spm.stats.fmri_spec.mthresh = 0.8;
matlabbatch{3}.spm.stats.fmri_spec.mask = {''};
matlabbatch{3}.spm.stats.fmri_spec.cvi = 'AR(1)';
%--------------------------------------------------------------------------
% GLM PPImodel估计
%--------------------------------------------------------------------------
matlabbatch{4}.spm.stats.fmri_est.spmmat(1) = cfg_dep('fMRI model specification: SPM.mat File', substruct('.','val', '{}',{3}, '.','val', '{}',{1}, '.','val', '{}',{1}), substruct('.','spmmat'));
matlabbatch{4}.spm.stats.fmri_est.write_residuals = 0;
matlabbatch{4}.spm.stats.fmri_est.method.Classical = 1;
%--------------------------------------------------------------------------
% PPI contrast
%--------------------------------------------------------------------------
matlabbatch{5}.spm.stats.con.spmmat(1) = cfg_dep('Model estimation: SPM.mat File', substruct('.','val', '{}',{4}, '.','val', '{}',{1}, '.','val', '{}',{1}), substruct('.','spmmat'));
matlabbatch{5}.spm.stats.con.consess{1}.tcon.name = 'PPI_interaction';
matlabbatch{5}.spm.stats.con.consess{1}.tcon.weights = [1 0 0 0];
matlabbatch{5}.spm.stats.con.consess{1}.tcon.sessrep = 'none';
matlabbatch{5}.spm.stats.con.delete = 0;
%--------------------------------------------------------------------------
% SPM run
%--------------------------------------------------------------------------
spm_jobman('run',matlabbatch);
end

五、批量移动结果,准备统计分析

%时间:20220919
%程序员:wishing
%批量复制结果文件
%--------------------------------------------------------------------------
% 设置循环
clear,clc;
%————————————————————————————————
for i=[2:4,6:7,11:13,15:22,24:31,34:37]
p(i) = string(num2str(i, "%02d"));%生成“01”“02”一直到“37”的字符串,剔除1, 5, 8, 9, 10, 14, 23, 32, 33)
end
p = rmmissing(p);%删除跳过产生的缺失值
%----------------------------------------------------------------
% 设置文件夹
%————————————————————————————————
for n=1:length(p)
    sourcepath=string(append('D:\desk\speed\data_1\Dp_spm\PPI\RSBA18\RW\sub_1',p(n),'\con_0001.nii'));
    targetpath=string(append('D:\desk\speed\data_1\Dp_spm\PPI\RSBA18\RW\RW28CON\sub_1',p(n),'.nii'));
    copyfile(sourcepath,targetpath);
end
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 206,602评论 6 481
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 88,442评论 2 382
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 152,878评论 0 344
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 55,306评论 1 279
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 64,330评论 5 373
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 49,071评论 1 285
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 38,382评论 3 400
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 37,006评论 0 259
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 43,512评论 1 300
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,965评论 2 325
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,094评论 1 333
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,732评论 4 323
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 39,283评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,286评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,512评论 1 262
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 45,536评论 2 354
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,828评论 2 345

推荐阅读更多精彩内容