GroupStatsSPM/jpa_getSigResults.m at develop · jpalbrecht/GroupStatsSPM · GitHub

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
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
function sigRes = jpa_getSigResults(TabDat, atlas, pValue, varargin)
% Function that filters significant results from a SPM-Matlab
% TabDat-Structure.
%
% Syntax:
%    jpa_getSigResults(TabDat, atlas, pValue [, cFWE, cFDR, pFWE, pFDR])
%
% Inputs:
%    TabDat  - SPM-Matlab Structure containing the T/F-Values of a
%               statistical Test
%    atlas   - Path to a file containing a anatomical region-name to every
%               Voxel-Coordinate
%    pValue  - p-Value
%    cFWE    - Boolean, enables the Cluster-Family-Wise-Error-Correction
%               evaluation
%    cFDR    - Boolean, enables the Cluster-False Detection Rate evaluation
%    pFWE    - Boolean, enables the Peak-Family-Wise-Error-Correction
%               evaluation
%    pFDR    - Boolean, enables the Peak-False Detection Rate evaluation
%
% Outputs:
%    sigRes    - Struct containing all Values of TabDat wich are over a
%                   certain pValue, labeled with anatomical region
%
% Example:
%           jpa_getSigResults(TabDat, 'C:/atlas1',0.05, 1, 1, 1, 1)
%
% Other m-files required: none
% Subfunctions: getSigRes
% MAT-files required: none
%
% See also:

% Author: Jan Albrecht
% Work address:
% email: jan-philipp.albrecht@charite.de, j.p.albrecht@fu-berlin.de
% Website:
% Sep 2015; Last revision: 29-Sep-2015

%------------- BEGIN CODE --------------

%% Input Parameters check
if nargin < 3
    error('Wrong number of Arguments');
elseif nargin == 3
    cFWE  = 1;
    cFDR = 1;
    pFWE = 1;
    pFDR = 1;
elseif nargin == 4
    cFWE  = varargin{1};
    cFDR = 1;
    pFWE = 1;
    pFDR = 1;
elseif nargin == 5
    cFWE  = varargin{1};
    cFDR = varargin{2};
    pFWE = 1;
    pFDR = 1;
elseif nargin == 6
    cFWE  = varargin{1};
    cFDR = varargin{2};
    pFWE = varargin{3};
    pFDR = 1;
else
    cFWE  = varargin{1};
    cFDR = varargin{2};
    pFWE = varargin{3};
    pFDR = varargin{4};
end

%% Initialization
% Initialisation Cluster FWE
sigRes.cFWE.dat = [];
sigRes.cFWE.peakCoord = [];
sigRes.cFWE.label = [];
% Initialisation Cluster FDR
sigRes.cFDR.dat = [];
sigRes.cFDR.peakCoord = [];
sigRes.cFDR.label = [];
% Initialisation Peak FWE
sigRes.pFWE.dat = [];
sigRes.pFWE.peakCoord = [];
sigRes.pFWE.label = [];
% Initialisation Peak FDR
sigRes.pFDR.dat = [];
sigRes.pFDR.peakCoord = [];
sigRes.pFDR.label = [];

%% run
% check for emptyness
if ~isfield(TabDat,'dat') || isempty(TabDat.dat)
    return;
end
try
    % try to load atlas
    xA = spm_atlas('load',atlas);
catch THEM
    disp('Could not load Atlas! because of reason:');
    disp(THEM.message);
    disp('set label to "NONE"');
    xA = 'NONE';
end

% get peakCoords
peakCoords = TabDat.dat(:,end);
% check sig Results for Cluster FWE
if cFWE
    % get FWE out of Table
    [sigRes.cFWE.dat, sigRes.cFWE.peakCoord, sigRes.cFWE.label] = getSigRes(TabDat,pValue,peakCoords,xA,3);
end
% check sig Results for Cluster FDR
if cFDR
    % get FDR out of Table
    [sigRes.cFDR.dat, sigRes.cFDR.peakCoord, sigRes.cFDR.label] = getSigRes(TabDat,pValue,peakCoords,xA,4);
end


% check sig Results for peak FWE
if pFWE
    [sigRes.pFWE.dat, sigRes.pFWE.peakCoord, sigRes.pFWE.label] = getSigRes(TabDat,pValue,peakCoords,xA,7);
end
% check sig Results for peak FDR
if pFDR
    [sigRes.pFDR.dat, sigRes.pFDR.peakCoord, sigRes.pFDR.label] = getSigRes(TabDat,pValue,peakCoords,xA,8);
end

end

% Function that filters the data to a certain threshold and gets the
% anatomical region of Coords
function [filter, peak, label] = getSigRes(TabDat,pValue,peakCoords,xA,TableNR)
raw = cell2mat(TabDat.dat(:,TableNR));
% filter
filter = raw(raw < pValue);
peak = peakCoords(raw < pValue);
% label
label = cell.empty;
for i=1:1:length(peak)
    if ~ischar(xA)
        label{i,1} = spm_atlas('query',xA,peak{i});
    else
        label{i,1} = 'NONE';
    end
end
end
%------------- END CODE --------------