Browse code

coordTabel works for JZ006

git-svn-id: https://svn.discofish.de/MATLAB/spmtoolbox/SVMCrossVal@138 83ab2cfd-5345-466c-8aeb-2b2739fb922d

Christoph Budziszewski authored on 25/02/2009 15:49:55
Showing 9 changed files
... ...
@@ -1,36 +1,4 @@
1
-% filenameList as in SPM.xY.VY
2
-% voxelList in [x y z], 1 coordinate per row
3
-
4 1
 function extr = calculateImageData(filenameList, voxelList)
5
-TRANSFORM_METHOD = 'single';
6
-%function extr = calculateImageData(voxelList,des,smoothed)
7
-% global USE_DRIVE_CHECK_HACK;
8
-
9
-% dtype='PSTH';
10
-
11
-% switch dtype 
12
-%     case 'PSTH'
13
-%         V=des.xY.VY;
14
-%     case 'betas'
15
-%         V=des.Vbeta;
16
-% end;
17
-
18
-% if USE_DRIVE_CHECK_HACK
19
-% if V(1).fname(1)~='D'
20
-%      for z=1:length(V) % Change Drive Letter - HACK!
21
-%       	V(z).fname(1) = 'D';
22
-%      end; 
23
-% end
24
-% end
25
-
26
-% if (~smoothed)
27
-%   for z=1:length(V) % Change smoothed Filename - HACK!
28
-%       % D:....SUBJECTID\session\swfanders...
29
-%       % D:....SUBJECTID\session\wfanders...
30
-%       tmp = findstr(filesep,V(z).fname);
31
-%       V(z).fname = strcat(V(z).fname(1:tmp(end)),V(z).fname(tmp(end)+2:end));
32
-%   end;
33
-% end
34 2
 
35 3
 V = filenameList;
36 4
 
... ...
@@ -39,66 +7,38 @@ nVoxel = size(vox,1);
39 7
 nImage = numel(V);
40 8
 
41 9
 for kImage=1:nImage
42
-% 	extr(kImage) = struct(...
43
-%         'dat',  ...
44
-%         'val',   repmat(NaN, [1 nVoxel]),...
45
-% 		'mean',  repmat(NaN, [1 nVoxel]),...
46
-% 		'sum',   repmat(NaN, [1 nVoxel]),...
47
-% 		'nvx',   repmat(NaN, [1 nVoxel]),...
48
-% 		'posmm', repmat(NaN, [3 nVoxel]),...
49
-% 		'posvx', repmat(NaN, [3 nVoxel]));
10
+    roicenter = round(inv(V(kImage).mat)*[vox, ones(nVoxel,1)]');
11
+    x = roicenter(1,:);
12
+    y = roicenter(2,:);
13
+    z = roicenter(3,:);
50 14
 
51
-    switch TRANSFORM_METHOD
52
-        case 'single'
53
-            roicenter = round(inv(V(kImage).mat)*[vox, ones(nVoxel,1)]');
54
-            x = roicenter(1,:);
55
-            y = roicenter(2,:);
56
-            z = roicenter(3,:);
57
-        case 'image'
58
-             x = []; y = []; z = [];
59
-            for iRoiFile = 1:nRoiFiles 
60
-                [x1 y1] = ndgrid(1:V(k).dim(1),1:V(k).dim(2));
61
-                for p = 1:V(k).dim(3) % resample mask Vm(iRoiFile) in space of V(k)
62
-                    B = spm_matrix([0 0 -p 0 0 0 1 1 1]);
63
-                    M = inv(B*inv(V(k).mat)*Vm(iRoiFile).mat);
64
-                    msk = find(spm_slice_vol(Vm(iRoiFile),M,V(k).dim(1:2),0));
65
-                    if ~isempty(msk)
66
-                        z1 = p*ones(size(msk(:)));
67
-                        x = [x; x1(msk(:))];
68
-                        y = [y; y1(msk(:))];
69
-                        z = [z; z1];
70
-                    end;
71
-                end;
72
-            end
73
-    end
74
-    
75 15
 
76
-	for iVoxel = 1:nVoxel
16
+    for iVoxel = 1:nVoxel
77 17
 
78
-%         if rad==0
18
+        if rad==0
79 19
             x = roicenter(1,iVoxel);
80 20
             y = roicenter(2,iVoxel);
81 21
             z = roicenter(3,iVoxel);
82
-%         else
83
-%             tmp = spm_imatrix(V(kImage).mat);
84
-%             vdim = tmp(7:9);
85
-%             vxrad = ceil((rad*ones(1,3))./(ones(nVoxel,1)*vdim))';
86
-%             [x y z] = ndgrid(-vxrad(1,iVoxel):sign(vdim(1)):vxrad(1,iVoxel), ...
87
-%                       -vxrad(2,iVoxel):sign(vdim(2)):vxrad(2,iVoxel), ...
88
-%                       -vxrad(3,iVoxel):sign(vdim(3)):vxrad(3,iVoxel));
89
-%             sel = (x./vxrad(1,iVoxel)).^2 + (y./vxrad(2,iVoxel)).^2 + ...
90
-%                   (z./vxrad(3,iVoxel)).^2 <= 1;
91
-%             x = roicenter(1,iVoxel)+x(sel(:));
92
-%             y = roicenter(2,iVoxel)+y(sel(:));
93
-%             z = roicenter(3,iVoxel)+z(sel(:));
94
-%         end;
95
-
96
-
97
-		dat = spm_sample_vol(V(kImage), x, y, z,0);
22
+        else
23
+            tmp = spm_imatrix(V(kImage).mat);
24
+            vdim = tmp(7:9);
25
+            vxrad = ceil((rad*ones(1,3))./(ones(nVoxel,1)*vdim))';
26
+            [x y z] = ndgrid(-vxrad(1,iVoxel):sign(vdim(1)):vxrad(1,iVoxel), ...
27
+                -vxrad(2,iVoxel):sign(vdim(2)):vxrad(2,iVoxel), ...
28
+                -vxrad(3,iVoxel):sign(vdim(3)):vxrad(3,iVoxel));
29
+            sel = (x./vxrad(1,iVoxel)).^2 + (y./vxrad(2,iVoxel)).^2 + ...
30
+                (z./vxrad(3,iVoxel)).^2 <= 1;
31
+            x = roicenter(1,iVoxel)+x(sel(:));
32
+            y = roicenter(2,iVoxel)+y(sel(:));
33
+            z = roicenter(3,iVoxel)+z(sel(:));
34
+        end;
35
+
36
+
37
+        dat = spm_sample_vol(V(kImage), x, y, z,0);
98 38
         extr(kImage).dat(iVoxel)      = dat;
99
-		extr(kImage).mean(iVoxel)     = nanmean(dat);
39
+        extr(kImage).mean(iVoxel)     = nanmean(dat);
100 40
         extr(kImage).nvx(iVoxel)      = numel(dat);
101
-	end;
41
+    end;
102 42
 
103 43
 end;
104 44
 end
105 45
\ No newline at end of file
106 46
new file mode 100644
... ...
@@ -0,0 +1,45 @@
1
+%% subject loop
2
+function decode = claculateMultiSubjectDecodePerformance(timelineParams,calculateParams,paramModel)
3
+
4
+decode = struct;
5
+decode.decodePerformance = [];
6
+decode.rawTimeCourse     = [];
7
+
8
+for subjectCell = calculateParams.subjectSelection
9
+    SubjectID = cell2mat(subjectCell);
10
+    namehelper = strcat('s',SubjectID); %Vars can not start with numbers.
11
+
12
+    display('loading SPM.mat ...');
13
+    spm = load(fullfile(calculateParams.PROJECT_BASE_PATH,SubjectID,calculateParams.PROJECT_RESULT_PATH));
14
+    display('... done.');
15
+
16
+    % per subject params..
17
+    subjectParams = struct;
18
+    
19
+    subjectParams.des             = spm.SPM;
20
+
21
+
22
+    subjectParams.voxelList       = mapVoxelList(voxelList,SubjectID);
23
+
24
+    
25
+    subjectParams.SubjectID       = SubjectID;
26
+    subjectParams.namehelper      = namehelper;
27
+%     assignin('base','calculateParams',calculateParams);
28
+
29
+    display(sprintf('calculating cross-validation performance time-shift for Subject %s. Please Wait. ...',SubjectID));
30
+    display('switching off all warnings');
31
+    warning_state               = warning('off','all');
32
+    display('calculating ...');
33
+    decode.(namehelper)         = calculateDecodePerformance(timelineParams,calculateParams,subjectParams);
34
+
35
+    display('... done');
36
+    display('restoring warnings');
37
+    warning(warning_state);
38
+
39
+    decode.decodePerformance    = [decode.decodePerformance decode.(namehelper).decodePerformance];
40
+    decode.rawTimeCourse        = [decode.rawTimeCourse decode.(namehelper).rawTimeCourse];
41
+
42
+    assignin('base','decode',decode);
43
+end
44
+
45
+end
0 46
new file mode 100644
... ...
@@ -0,0 +1,19 @@
1
+function scoords = getSubjectCoordinates(gcoords,map)
2
+m = java.util.HashMap ;
3
+
4
+for i = 1:numel(map.names)
5
+    key = cell2mat(map.names(i));
6
+    val = map.coords(i,:);
7
+    m.put(key,val);
8
+end
9
+
10
+scoords = [];
11
+
12
+for c = 1:numel(gcoords)
13
+    key = strtrim(gcoords(c).name);
14
+    val = m.get(key)';
15
+    mod = gcoords(c).mod;
16
+    scoords = [scoords; val+mod];
17
+end
18
+
19
+end
0 20
\ No newline at end of file
... ...
@@ -5,23 +5,35 @@ disp('RUN');
5 5
 timeLine = getTimeLineParams(model);
6 6
 subjects = getSubjectCellList(model);
7 7
 classDef = parseClassDef(model);
8
+
9
+mask     = ['^' cell2mat(getImageFileMask(model)) '.*\.img$'];
8 10
 % images
9 11
 % normalization
10 12
 
11 13
 switch task
12 14
     case 'COORD'
13 15
         disp('COORD');
14
-        coordinates= 'parse me';
15
-        runCoordTable()
16
+        coordargs = struct;
17
+        coordargs.subjects = subjects;
18
+        coordargs.timeline = timeLine;
19
+        coordargs.basedir = model.baseDir;
20
+        coordargs.sessionList = 1:3;
21
+        coordargs.eventList = classDef.eventMatrix;
22
+        coordargs.coords = parseCoordinateTextField(model);
23
+        coordargs.mask   = mask;
24
+        
25
+        runCoordTable(coordargs)
16 26
         
17 27
     case 'ROI'
18 28
         disp('ROI');
19 29
         roiargs = struct;
20
-        roiargs.subjects = subjects;
21
-        roiargs.timeline = timeLine;
22
-        roiargs.classes  = classDef;
23
-        roiargs.mask     = getImageFileMask(model);
24
-        roiargs.basedir  = model.baseDir;
30
+        roiargs.subjects    = subjects;
31
+        roiargs.timeline    = timeLine;
32
+        roiargs.classes     = classDef;
33
+        roiargs.mask        = mask;
34
+        roiargs.basedir     = model.baseDir;
35
+        roiargs.sessionList = 1:3;
36
+        roiargs.eventList   = classDef.eventMatrix;
25 37
         
26 38
         
27 39
         assignin('base','roiargs',roiargs);
... ...
@@ -1,5 +1,5 @@
1
-function voxelList = parseCoordinateTextField(txtFieldHandel)
2
-s = get(txtFieldHandel,'String');
1
+function voxelList = parseCoordinateTextField(model)
2
+s = get(model.txtVoxelDef,'String');
3 3
 rows = size(s,1);
4 4
 voxelList = [];
5 5
     for i = 1:rows 
6 6
deleted file mode 100644
... ...
@@ -1,51 +0,0 @@
1
-function voxelList = parseVoxelList(model,multisubjectid)
2
-        voxelList = [];
3
-
4
-        %<ROI Name>,<ROI Modifier>;
5
-        txt = get(model.txtVoxelDef,'String');
6
-        map = model.subjectMap;
7
-
8
-        switch nargin
9
-            case 1
10
-                SubjectID = getSubjectIDString(model);
11
-            case 2
12
-                SubjectID = multisubjectid;
13
-            otherwise
14
-                error('spmtoolbox:SVMCrossVal:parseVoxelList:nargin','wrong number of arguments given');
15
-        end
16
-         
17
-        rows  = size(txt,1);
18
-        
19
-        for i = 1:rows 
20
-            if all(isspace(txt(i,:)))
21
-                continue;
22
-            end
23
-            line = txt(i,:);
24
-            roi = parseROIName(line);
25
-            roimod = parseModifier(line);
26
-            voxelList = [voxelList; getCoordinate(map,SubjectID,roi)+eval(roimod)];
27
-        end
28
-end
29
-
30
-function roi = parseROIName(line)
31
-c = 1;
32
-roi = '';
33
-while line(c)~=','
34
-    roi = [roi line(c)];
35
-    c = c+1;
36
-end
37
-roi = strtrunc(roi);
38
-end
39
-
40
-function roimod = parseModifier(line)
41
-c = 1;
42
-while line(c)~=',' % skip roi
43
-    c=c+1;
44
-end
45
-c=c+1; % skip ','
46
-roimod='';
47
-while line(c) ~= ';'
48
-    roimod = [roimod line(c)];
49
-    c = c+1;
50
-end
51
-end
... ...
@@ -1,98 +1,53 @@
1
-function ppsth = runCoordTable(baseDir,subjects,timeline,coordinates,volumeFileList,classDef)
2
-%         paramModel = varargin{1};
1
+function runCoordTable(args)
3 2
     disp('run coord table')
4 3
     
5
-    disp('does not work yet.')
6
-    return
7
-        
8
-        PROJECT_BASE_PATH = baseDir;
9
-        PROJECT_RESULT_PATH = 'results\SPM.mat';
10
-
4
+    subjects = args.subjects;
5
+    nSubjects = size(subjects);
6
+    sessionlist = args.sessionList;
11 7
 
12
-        timelineParams = timeline;
13
-        
14
-        % common params
15
-        calculateParams  = struct;
16
-        
17
-        calculateParams.PROJECT_BASE_PATH    = PROJECT_BASE_PATH;
18
-        calculateParams.PROJECT_RESULT_PATH  = PROJECT_RESULT_PATH;
19
-        
20
-        calculateParams.RANDOMIZE       = 0;
21
-        calculateParams.sessionList     = 1:3;
22
-        
23
-        calculateParams.eventList       = classStruct.eventMatrix; %[9,11,13; 10,12,14];
8
+    for s = 1:nSubjects
9
+        subjectStruct{s}.dir = fullfile(args.basedir,cell2mat(subjects(s)));
10
+        d = load(fullfile(subjectStruct{s}.dir,'results','SPM.mat'));
11
+        subjectStruct{s}.des = d.SPM;
12
+        subjectStruct{s}.name = cell2mat(subjects(s));
24 13
         
25
-        subjectSelection = subjects;
26
-        calculateParams.subjectSelection = subjectSelection;
14
+        map = load(fullfile(subjectStruct{s}.dir,'results','roi','coord_map.mat'));
15
+        subjectStruct{s}.coords = getSubjectCoordinates(args.coords,map);
27 16
         
28
-        decode = claculateMultiSubjectDecodePerformance(timelineParams,calculateParams,paramModel);
17
+        disp('fetching volume definitions, please wait');
18
+        subjectStruct{s}.volumes = spm_vol(getImageFileList(subjectStruct{s}.dir,sessionlist,args.mask));
29 19
 
30
-        display('Finished calculations.');
31
-        display('Plotting...');
32
-
33
-        plotParams                   = struct;
34 20
         
35
-%         plotParams.SVMCROSSVAL_CROSSVAL_METHOD_DEF = SVMCROSSVAL_CROSSVAL_METHOD_DEF;
36
-        plotParams.CROSSVAL_METHOD     = calculateParams.CROSSVAL_METHOD;
21
+        rawData = calculateImageData(subjectStruct{s}.volumes,subjectStruct{s}.coords)
37 22
         
38
-        plotParams.nClasses          = length(calculateParams.classList);
39
-
40
-        plotParams.decodePerformance = decode.decodePerformance;
41
-        plotParams.rawTimeCourse     = decode.rawTimeCourse;
42
-        plotParams.SubjectID         = subjectSelection;
43
-        plotParams.smoothed          = boolToYesNoString(calculateParams.smoothed);
44
-
45
-        assignin('base','plotParams',plotParams);
46
-%         plotDecodePerformance(params.psthStart,params.psthEnd,params.nClasses,decode.decodeTable,params.frameShiftStart,params.frameShiftEnd,decode.rawTimeCourse);
47
-        plotDecodePerformance(timelineParams,plotParams);
48
-            
49
-        display('all done.');
23
+        disp('done');
24
+    end
50 25
 end
51
-    
52
-
53
-%% subject loop
54
-function decode = claculateMultiSubjectDecodePerformance(timelineParams,calculateParams,paramModel)
55
-
56
-decode = struct;
57
-decode.decodePerformance = [];
58
-decode.rawTimeCourse     = [];
59 26
 
60
-for subjectCell = calculateParams.subjectSelection
61
-    SubjectID = cell2mat(subjectCell);
62
-    namehelper = strcat('s',SubjectID); %Vars can not start with numbers.
63 27
 
64
-    display('loading SPM.mat ...');
65
-    spm = load(fullfile(calculateParams.PROJECT_BASE_PATH,SubjectID,calculateParams.PROJECT_RESULT_PATH));
66
-    display('... done.');
67 28
 
68
-    % per subject params..
69
-    subjectParams = struct;
70
-    
71
-    subjectParams.des             = spm.SPM;
72
-
73
-
74
-    subjectParams.voxelList       = mapVoxelList(voxelList,SubjectID);
75
-
76
-    
77
-    subjectParams.SubjectID       = SubjectID;
78
-    subjectParams.namehelper      = namehelper;
79
-%     assignin('base','calculateParams',calculateParams);
80
-
81
-    display(sprintf('calculating cross-validation performance time-shift for Subject %s. Please Wait. ...',SubjectID));
82
-    display('switching off all warnings');
83
-    warning_state               = warning('off','all');
84
-    display('calculating ...');
85
-    decode.(namehelper)         = calculateDecodePerformance(timelineParams,calculateParams,subjectParams);
86
-
87
-    display('... done');
88
-    display('restoring warnings');
89
-    warning(warning_state);
90
-
91
-    decode.decodePerformance    = [decode.decodePerformance decode.(namehelper).decodePerformance];
92
-    decode.rawTimeCourse        = [decode.rawTimeCourse decode.(namehelper).rawTimeCourse];
93
-
94
-    assignin('base','decode',decode);
95
-end
96
-
97
-end
98 29
 
30
+% 
31
+% %         decode = claculateMultiSubjectDecodePerformance(timelineParams,calculateParams,paramModel);
32
+% 
33
+%         display('Finished calculations.');
34
+%         display('Plotting...');
35
+% 
36
+%         plotParams                   = struct;
37
+%         
38
+% %         plotParams.SVMCROSSVAL_CROSSVAL_METHOD_DEF = SVMCROSSVAL_CROSSVAL_METHOD_DEF;
39
+%         plotParams.CROSSVAL_METHOD     = calculateParams.CROSSVAL_METHOD;
40
+%         
41
+%         plotParams.nClasses          = length(calculateParams.classList);
42
+% 
43
+%         plotParams.decodePerformance = decode.decodePerformance;
44
+%         plotParams.rawTimeCourse     = decode.rawTimeCourse;
45
+%         plotParams.SubjectID         = subjectSelection;
46
+%         plotParams.smoothed          = boolToYesNoString(calculateParams.smoothed);
47
+% 
48
+%         assignin('base','plotParams',plotParams);
49
+% %         plotDecodePerformance(params.psthStart,params.psthEnd,params.nClasses,decode.decodeTable,params.frameShiftStart,params.frameShiftEnd,decode.rawTimeCourse);
50
+%         plotDecodePerformance(timelineParams,plotParams);
51
+%             
52
+%         display('all done.');
53
+% 
... ...
@@ -1,14 +1,15 @@
1 1
 function runROIImageMaskMode(args)
2 2
 
3 3
 subjects = args.subjects;
4
-mask     = ['^' cell2mat(args.mask) '.*\.img$'];
5
-    
4
+   
6 5
 nSubjects = size(subjects);
7
-sessionlist = 1:3;
6
+sessionlist = args.sessionList;
8 7
 
9 8
 
10 9
 for s = 1:nSubjects
11 10
     subjectStruct{s}.dir = fullfile(args.basedir,cell2mat(subjects(s)));
11
+    d = load(fullfile(subjectStruct{s}.dir,'results','SPM.mat'));
12
+    subjectStruct{s}.des = d.SPM;
12 13
     subjectStruct{s}.name = cell2mat(subjects(s));
13 14
     subjectStruct{s}.roiFile = ui_selectRoiImage(...
14 15
         sprintf('Select ROI Files for %s',subjectStruct{s}.name),...
... ...
@@ -21,12 +22,22 @@ for s = 1:nSubjects
21 22
     % load image data
22 23
    
23 24
     disp('fetching volume definitions, please wait');
24
-    subjectStruct{s}.volumes = spm_vol(getImageFileList(subjectStruct{s}.dir,sessionlist,mask));
25
+    subjectStruct{s}.volumes = spm_vol(getImageFileList(subjectStruct{s}.dir,sessionlist,args.mask));
25 26
 
26 27
     disp('computing volume values, please wait');
27
-    subjectStruct{s}.rawData = calculateRoiImageData(subjectStruct{s}.volumes,subjectStruct{s}.roiFile);
28
+    
29
+    rawData = calculateRoiImageData(subjectStruct{s}.volumes,subjectStruct{s}.roiFile);
28 30
     % calculate psth
29
-          
31
+    pstopts.des = subjectStruct{s}.des;
32
+    pstopts.eventList = args.eventList;
33
+    pstopts.sessionList = sessionlist;
34
+    
35
+    assignin('base','pstopts',pstopts);
36
+    assignin('base','rawData',rawData);
37
+
38
+    disp('computing psth');
39
+    subjectStruct{s}.psth    = calculatePST(args.timeline,pstopts,rawData);
40
+    disp('done');
30 41
 end
31 42
 
32 43
 assignin('base','subjectStruct',subjectStruct);
... ...
@@ -281,18 +281,18 @@ function model = createFirstStepPanel(model,parent,DEFAULT)
281 281
         pNorm = uipanel(parent,'Title','Normalization','Position',cell2mat(main_grid(2,3)));
282 282
         set(pNorm,'BackgroundColor','w');
283 283
         
284
-        createLabel(pNorm,[0 0.75 1 0.25],'PST Normalization');
285
-        norm1Model = {'DUMMY norm A','DUMMY norm B'};
284
+        createLabel(pNorm,[0 0.75 1 0.25],'psth norm4SVM');
285
+        norm1Model = {'none','mean','minmax'};
286 286
         model.selNormPST = uicontrol(pNorm,'Style','popupmenu',...
287 287
             'Units','normalized',...
288 288
             'Position',[0.0 0.5 1 0.25],...
289 289
             'String',norm1Model,...
290 290
             'UserData',norm1Model,...
291
-            'Value',1);
291
+            'Value',2);
292 292
          set(model.selNormPST,'BackgroundColor','w');   
293 293
         
294
-        createLabel(pNorm,[0 0.25 1 0.25],'Class-Grouping Normalization');
295
-        norm2Model = {'DUMMY norm X','DUMMY norm Y'};
294
+        createLabel(pNorm,[0 0.25 1 0.25],'Col Bias removal');
295
+        norm2Model = {'on','off'};
296 296
         model.selNormClass = uicontrol(pNorm,'Style','popupmenu',...
297 297
             'Units','normalized',...
298 298
             'Position',[0.0 0.0 1 0.25],...
... ...
@@ -313,7 +313,7 @@ function model = createFirstStepPanel(model,parent,DEFAULT)
313 313
         btnRunButton2 = uicontrol(pButtons,'String','run full Brain Searchlight',...
314 314
             'Units','normalized','Position',[0.33 0 0.33 1]);
315 315
         set(btnRunButton2,'Callback',{@cbRunPreprocessing,model,'FBS'}); % set here, because of model.  
316
-        set(btnRunButton2,'Enable','on');
316
+        set(btnRunButton2,'Enable','off');
317 317
         
318 318
         btnRunButton3 = uicontrol(pButtons,'String','run ROI-Image processing',...
319 319
             'Units','normalized','Position',[0.66 0 0.33 1]);