-
Notifications
You must be signed in to change notification settings - Fork 0
/
ASLOutlierRemoval.m
73 lines (59 loc) · 5.1 KB
/
ASLOutlierRemoval.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
function SUBJECT = ASLOutlierRemoval(SUBJECT, prefix, method)
% ClinicalASL toolbox 2023, JCWSiero
tic
outputmap = [SUBJECT.ASLdir prefix '_BASIL_perDYNAMIC/'];
if ~exist(outputmap,'dir')
mkdir(outputmap)
end
DIMS = size(SUBJECT.(prefix).ASL_label1label2_allPLD);
if ~exist([outputmap 'DYNAMIC_1'],'dir')
% Loop: compute CBF per dynamics, making sure we have 2 (control,label) x nPLD sized datasets
for i=1:SUBJECT.NREPEATS
%save data for BASIL analysis
SaveDataNII(reshape(SUBJECT.(prefix).ASL_label1label2_allPLD(:,:,:,2*i-1:2*i,:),DIMS(1),DIMS(2),DIMS(3), SUBJECT.NPLDS*2), [outputmap prefix '_oneDYNAMIC_label1label2'], SUBJECT.dummyfilenameSaveNII, 1, [],SUBJECT.TR)
% perform BASIL analysis
locationASL_multiPLD = [outputmap prefix '_oneDYNAMIC_label1label2'];
locationM0 = [SUBJECT.ASLdir prefix '_M0'];
locationMask = [SUBJECT.ASLdir prefix '_M0_brain_mask'];
disp(['**************************************************************************************** BASIL analysis on DYNAMIC = ' num2str(i) ' *****************************************************'])
ASLBASILanalysis(SUBJECT, locationASL_multiPLD, locationM0, locationMask, [outputmap 'DYNAMIC_' num2str(i)], [1:SUBJECT.NPLDS], SUBJECT.locationBASILinfo, 'artoff')
SUBJECT.(prefix).CBF_DYNAMIC(:,:,:,i) = double(niftiread([SUBJECT.ASLdir prefix '_BASIL_perDYNAMIC/DYNAMIC_' num2str(i) '/native_space/perfusion_calib']));
end
% save CBF per dynamic
SaveDataNII(SUBJECT.(prefix).CBF_DYNAMIC,[SUBJECT.ASLdir prefix '_CBF_perDYNAMIC'], SUBJECT.dummyfilenameSaveNII ,1,[],SUBJECT.TR);
disp(['Finished saving CBF data per dynamic: ' prefix]);
else
disp(['Seems BASIL analysis per Dynamic has already been done, going straight to performing ASL oulierremoval by Duloi et al...' newline])
disp(['Loading previous BASIL analysis per dynamic...'])
for i=1:SUBJECT.NREPEATS
SUBJECT.(prefix).CBF_DYNAMIC(:,:,:,i) = double(niftiread([SUBJECT.ASLdir prefix '_BASIL_perDYNAMIC/DYNAMIC_' num2str(i) '/native_space/perfusion_calib']));
end
end
% %% ASL OutlierRemoval (volumes), SCORE method by Duloi et al JMRI 2017 %%%
[IOR_allsteps, IOR_step1, IOR_step2, NoOR_logical] = ASLOutlierRemovalPerform(SUBJECT.(prefix).CBF_DYNAMIC, SUBJECT.(prefix).brainmask, SUBJECT.outlierFactor, SUBJECT.(prefix).GMmask, SUBJECT.(prefix).WMmask, SUBJECT.(prefix).CSFmask, method); %method: 'OnlyHighCBF': only use step1 of outlier removal (high CBF volumes), 'Duloi': or step1 + step2 (Duloi)
if isempty(IOR_allsteps)
IOR_allsteps=0;
end
% Write identified outliers to .txt
writematrix([numel(IOR_step1), numel(IOR_step2)],[SUBJECT.ASLdir prefix '_OutlierInformation_nvols_step1step2.txt'],'delimiter',' ');
writematrix(IOR_allsteps,[SUBJECT.ASLdir prefix '_OutlierInformation_whichvols_step1step2.txt'],'delimiter',' ');
%%%%%%%%%%%%%%%%%%%%%%%%%% Outlier removal and save to SUBJECT struct %%%%%%%%%%%%%%%%%%%%%%
SUBJECT.(prefix).NoOutliers_logical = NoOR_logical;
SUBJECT.(prefix).Ioutlier_allsteps = IOR_allsteps;
SUBJECT.(prefix).NREPEATS_OR = sum(SUBJECT.(prefix).NoOutliers_logical);
SUBJECT.(prefix).CBF_DYNAMIC_OR = SUBJECT.(prefix).CBF_DYNAMIC(:,:,:,SUBJECT.(prefix).NoOutliers_logical);
% save CBF outlier removed per dynamic
SaveDataNII(SUBJECT.(prefix).CBF_DYNAMIC_OR,[SUBJECT.ASLdir prefix '_CBF_perDYNAMIC_OR'], SUBJECT.dummyfilenameSaveNII ,1,[],SUBJECT.TR);
% remove outliers
NoOutliers_logical_forinterleaved = logical(zeros(1,2*SUBJECT.NREPEATS));
NoOutliers_logical_forinterleaved(1:2:end)=SUBJECT.(prefix).NoOutliers_logical; % odd indices
NoOutliers_logical_forinterleaved(2:2:end)=SUBJECT.(prefix).NoOutliers_logical; % even indices
SUBJECT.(prefix).ASL_label1label2_allPLD_OR = SUBJECT.(prefix).ASL_label1label2_allPLD(:,:,:,NoOutliers_logical_forinterleaved,:);
disp('Saving ASL data interleaved label control - OUTLIER REMOVED: all PLDs')
SaveDataNII(reshape(SUBJECT.(prefix).ASL_label1label2_allPLD_OR, DIMS(1),DIMS(2), DIMS(3), SUBJECT.NPLDS*SUBJECT.(prefix).NREPEATS_OR*2), [SUBJECT.ASLdir prefix '_allPLD_label1label2_OR'], SUBJECT.dummyfilenameSaveNII, 1,[], SUBJECT.TR) % save interleaved control label and for all PLDs as 4th dimension
disp('Saving ASL data interleaved label control - OUTLIER REMOVED: 2-to-last PLDs for CBF maps, "free" of arterial transit artefacts')
SaveDataNII(reshape(SUBJECT.(prefix).ASL_label1label2_allPLD_OR(:,:,:,:,2:end), DIMS(1),DIMS(2), DIMS(3), (SUBJECT.NPLDS-1)*SUBJECT.(prefix).NREPEATS_OR*2), [SUBJECT.ASLdir prefix '_2tolastPLD_label1label2_OR'], SUBJECT.dummyfilenameSaveNII, 1,[], SUBJECT.TR) % save interleaved control label and for all PLDs as 4th dimension
disp('Saving ASL data interleaved label control - OUTLIER REMOVED: 1-to-2 PLDs for ATA, aCBF maps, and spatial COV ')
SaveDataNII(reshape(SUBJECT.(prefix).ASL_label1label2_allPLD_OR(:,:,:,:,1:2), DIMS(1),DIMS(2), DIMS(3), 2*SUBJECT.(prefix).NREPEATS_OR*2), [SUBJECT.ASLdir prefix '_1to2PLD_label1label2_OR'], SUBJECT.dummyfilenameSaveNII, 1,[], SUBJECT.TR) % save interleaved control label and for all PLDs as 4th dimension
toc
end