MATLAB script for correlation maps

From Psych 221 Image Systems Engineering
Revision as of 09:19, 9 December 2009 by imported>Psych 204 (Created page with '% want to make a parameter map on subjectB's brain that shows the % correlation in each voxel for a particular voxel chosen from subjectA's % brain % load some data: clear all f…')
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search

% want to make a parameter map on subjectB's brain that shows the % correlation in each voxel for a particular voxel chosen from subjectA's % brain

% load some data: clear all format long g


%set directory for the two subjects to be compared % make sure you run this file from cjimmy if ispc

   subjADir='Y:\New_Localizers\dr021108-1p5mm-3mm';
   subjBDir='Y:\New_Localizers\kw011708-1p5mm-3mm';

else

   %     kevin subjects
   %     subjADir='/biac2/kgs/projects/New_Localizers/dr021108-1p5mm-3mm';
   % %     subjBDir='/biac2/kgs/projects/New_Localizers/kw011708-1p5mm-3mm';
   %     subjBDir='/biac2/kgs/projects/New_Localizers/kgs020408-1p5mm-3mm';
   % location of code
   home = '/biac2/kgs/projects/SameDifferent/'
   % nathan's prosopagnosic subjects

% subjADir='/biac2/kgs/projects/Prosopagnosia/fMRIData/112909whFMRI'; % subjADir='/biac2/kgs/projects/Prosopagnosia/fMRIData/020209_mr_fmri'; % subjBDir='/biac2/kgs/projects/Prosopagnosia/fMRIData/52309AF'; % subjADir = % '/biac2/kgs/projects/Kids/fmri/localizer/adult_jc_27yo_052408';

subjADir = '/biac2/kgs/projects/Kids/fmri/localizer/adult_al_22yo_051108/';

% subjBDir = % '/biac2/kgs/projects/Kids/fmri/localizer/adult_kw_25yo_090308'; subjBDir = '/biac2/kgs/projects/Kids/fmri/localizer/adult_dy_25yo_041908';

end

% %choose an ROI in subject A and B % waldemar % roiA='lPPA'; % roiA='rFFA'; % alana % roiB='AllRois_Combined';


% janice % roiA='lPPA_IOvAC_p3'; roiA='lFFA_MBvAC_p3'; % kevin roiB = 'allROIS';


%want to get the motion corrected data type which is usually 3 but might be %anything depending on the experimenter and experiments. dt = 3; %'MotionComp_RefScan5'; %then pick the scans we will use. % mtloc,eccbias1,eccbias2,4xretinotopy,objloc1,objloc2 % so to get our scans in order we have % scan = [8 9 1 2 3]; scan = [1 2 3 4];


% .. for subject A brain %go to subjects data directory cd(subjADir); %initialize our hidden inplane set to the right data type,scan(s), and roi % hiA = initHiddenInplane(dt, scan, roisA); hiA = initHiddenInplane(dt, scan);


% get the roi data from subject 1

%then we grab the multivoxel data from the ROI %this should give us mvA which is a struct with lots of useful stuff %in particular %mvA.trials has information about conditions and order of stimuli etc %mvA.coord has the coordinates for each voxel in the roi should be 3 x %numVoxels %mvA.tSeries has the time series for each voxel in % signal change %mvA.roi had roi info mvA = mv_init(hiA, roiA, scan, dt,1);

%then we apply our glm to each of the voxels in the roi %we will probably want to set these parameters ourselves so that we know %what happened %this gives us %mvA.glm has the betas, the design matrix and so on %for a block design there will be as many betas as there are conditions %for an event related design there will be timepoints x conditions x voxels %betas % mvA{i} = mv_applyGlm(mvA{i});


% now do the same for subject B % set scan numbers for B % scan=[1 2 3 4 5]; scan=[1 2 3 4];

% .. for all ROIs in subj B %instead invoke hidden inplane with multiple rois. then loop through the %rois when comparing to subjectA %go to subject b directory cd(subjBDir); % get dataTpe for subjectB dt=3; %load data with all rois hiB = initHiddenInplane(dt,scan,roiB); mvB = mv_init(hiB, roiB, scan, dt, 1); % mvB{i} = mv_applyGlm(mvB{i});



% can get the size of the brain map in the view for a particular scan by % doing

sz=dataSize(hiB,1); % which should be x,y, slice

% when we use mv_init on our rois we end up with coordinates in the field % mvA.roi.coords that have a particular numbering convention. we want to % transform that into the space that the view uses which we hope is done % with the function

% so we pick an arbitrary voxel from subject A's roi voxelInA = 10;

% then we compute the correlation with each voxel in subject B's roi


% vector for the correlations corrs=zeros(size(mvA.tSeries,2),size(mvB.tSeries,2)); coordinates=zeros(3,size(mvB.roi.coords,1)); for j = 1:size(mvA.tSeries, 2)

   for k = 1:size(mvB.tSeries, 2); %number of voxels in the ROI in B
       %get correlation between the time series of the two voxels
       r=corrcoef(mvA.tSeries(:,j), mvB.tSeries(:,k)); %records the correlations in a
       % record correlations
       corrs(j,k) = r(2);
       %     in order to get the same number of timeseries as anatomical
       %     coordinates must set preserveCoord variable in mv_init to be 1
       % some time series will be duplicated but we don't care
       %     coordinates(:,k) = mvB.roi.coords(:,k)
   end
   %

end

% get an average correlation

meancorr = mean(corrs,1);

% function [roiInd, coords] = roiIndices(view, coords, preserveCoords);


% in this case we think it should be [roiInd, coords] = roiIndices(hiB, mvB.roi.coords, 1);

% [oldcoords viewcoords] = roiIndices(view, ourroi, 0) % with the coords variable containing the mapping to the inplane.

% then we propose to make an array of the same size as the view ourmap = zeros(sz);

% then we will assign our correlations for each voxel to its location in % the matrix for i = 1:size(coords,2)

   ourmap(coords(1,i),coords(2,i),coords(3,i))=meancorr(i);

end

map = {ourmap,[]}; mapName = 'cor_meanallFFA_KWall'; mapUnits = 'r'; co ={[],[]};

pth = ['Inplane/GLMs/' mapName];

save(pth, 'map', 'mapName', 'mapUnits');

% save our parameter map


% -> load parameter map in mrvista and see what comes out.