Smoothing Code

%%Source: %%http://sites.google.com/site/mvlombardo/matlab-tutorials %%comments that don't start with 'Grace:' were written by website author

SubjArray={'F-1', 'F-2', 'F-4', 'I-6', 'H-5', 'C-1', 'C-4', 'C-5', 'C-7', 'J-5'}; %Grace: change subject names here! SubjArray2={'F1', 'F2', 'F4', 'I6', 'H5', 'C1', 'C4', 'C5', 'C7', 'J5'}; %Grace: ...and here, if your segmentation file's subject names were different Datamat=zeros(157, 189, 136, length(SubjArray));  % Grace: creates matrix where last number is subject index number

smoothmat_=zeros(157,189,136,length(SubjArray));

for i=1:length(SubjArray)

%datadir = ['/Users/grace/Desktop/PSYCH204A/' SubjArray{i} '/']; datadir = '/Users/grace/Desktop/PSYCH204A/'; fname = ['w' SubjArray2{i} '_NAcc_seg_allone.nii']; Data2Read = fullfile(datadir,fname);

% What I've done in Step 1 is to construct a filename such as /data/fMRI/12001.nii from its filepath (the variable called datadir) and filename (the variable called fname). Both these variables are string variables and I've used the function fullfile to essentially concatenate these two strings together. The end result is a string variable called Data2Read, which is essential a string that looks like this: /data/fMRI/12001.nii

% Step 2. Now pass Data2Read into the function spm_vol in order to read in the header information for the file. Note that the spm_vol function comes with the SPM toolbox in Matlab that is a standard software package for neuroimagers. Having SPM is essential for this part. HeaderInfo = spm_vol(Data2Read);

% In Step 2, the header information for the data you want to read in is now in your workspace. This is required for Step 3, where you will use the header information to actually read in the entire data.

% Step 3. Now use spm_read_vols to read in the data Data = spm_read_vols(HeaderInfo);

%% Grace: converting non zero matrix values to 1

indx=find(Data); Data(indx)=1; Datamat(:,:,:,i)=Data; %Grace: Datamat now contains all your matrix data

%% smoothing with smooth3

%smoothmat(:,:,:,i)=smooth3(Datamat(:,:,:,i), 'gaussian'); %smoothmat2(:,:,:,i)=smooth3(Datamat(:,:,:,i), 'gaussian', [1 1 1]); % smoothing kernel of 1, since voxel size is 1mm cube smoothmat3(:,:,:,i)=smooth3(Datamat(:,:,:,i), 'gaussian', [5 5 5]); %

for x=1:157

   for y=1:189
       for z=1:136
           if smoothmat3(x,y,z,i)>=0.3 % vary threshold here
               smoothmat_(x,y,z,i)=1;
           else
              smoothmat_(x,y,z,i)=0;
           end
       end
   end

end

%% % % % %% saving all values in file to '1' % % HeaderInfo.fname = ['w' SubjArray2{i} '_NAcc_seg_allten.nii'];  % This is where you fill in the new filename % % HeaderInfo.private.dat.fname = HeaderInfo.fname;  % This just replaces the old filename in another location within the header. % % % % % % Step 2. Now use spm_write_vol to write out the new data. You need to give spm_write_vol the new header information and corresponding data matrix % % spm_write_vol(HeaderInfo,Data); % end

combine_smooth_rethresh2=sum(smoothmat_,4); % % combinemat_g=sum(Datamat,4); %Grace: sum across subjects... now values indicate overlap across your subjects

%%

%%% %% histogram

index2=find(combine_smooth_rethresh2); %G: find non zero elements, so histogram doesn't have massive bar at 0; toplot=combine_smooth_rethresh2(index2); hist(toplot,9); % G: hardcoded to 9 bins, because maximum overlap was 9 ylabel('# voxels with x subjects overlapping') xlabel('x - number of subjects per voxel')

figure(2) rhist(toplot,9);

ylabel('relative frequency of voxels with x subjects overlapping') xlabel('x - number of subjects per voxel')

% %%% % % % % %% writing data out to nifti % % % % %%%% % % Step 1. Take the header information from a previous file with similar dimensions and voxel sizes and change the filename in the header. % HeaderInfo.fname = 'Combinesmooth_rethresh111_03.nii';  % This is where you fill in the new filename % HeaderInfo.private.dat.fname = HeaderInfo.fname;  % This just replaces the old filename in another location within the header. % % % Step 2. Now use spm_write_vol to write out the new data. You need to give spm_write_vol the new header information and corresponding data matrix % spm_write_vol(HeaderInfo,combine_smooth_rethresh2);  % where HeaderInfo is your header information for the new file, and Data is the image matrix corresponding to the image you'll be writing out.