MATLAB script for correlation maps
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 % roiA='lPPA'; % roiA='rFFA'; % alana % roiB='AllRois_Combined';
% roiA='lPPA_IOvAC_p3'; roiA='lFFA_MBvAC_p3';
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.