More Code

From Psych 221 Image Systems Engineering
Jump to navigation Jump to search

%% Analysis

% rgb color code for labels in ITK-Snap (match bar colors with labels)


% 0 0 0 0 0 0 0 "Clear Label" % 1 85 36 255 1 1 1 "Label 1" % 2 7 81 255 1 1 1 "Label 2" % 3 35 252 255 1 1 1 "Label 3" % 4 0 255 145 1 1 1 "Label 4" % 5 87 255 21 1 1 1 "Label 5" % 6 224 255 10 1 1 1 "Label 6" % 7 231 178 0 1 1 1 "Label 7" % 8 228 96 0 1 1 1 "Label 8" % 9 255 2 0 1 1 1 "Label 9" % 10 252 255 249 1 1 1 "Label 10" %

% Specify data sets with regular voxel size and double voxel size

dat = reg_vox_10;  % 1x1x1 voxel size, matrix for each subject dat_sum = sum(reg_vox_10,4);  % 1x1x1 voxel size, summed matrix across subjects

dat2 = MNIeach_ten;  % 2x2x2 voxel size, matrix for each subject dat2_sum = sum(MNIeach_ten,4);  % 2x2x2 voxel size, summed matrix across subjects dat_rep = MNI_NAcc_rep;  % SPM's representative subject NAc mask

%% 1x1x1 voxels

% Size - finds the # of voxels in each subjects mask then calculates the average and the s.e.

for k = 1:size(dat,4);  % 2x2x2

   indx = find(dat(:,:,:,k));
   subj_vox_count(k) = length(indx);
   

end

mean_nvox = mean(subj_vox_count); se_nvox = std(subj_vox_count)/(size(dat,4)^.5);


%% PLOTS

forhist = [];

for i=1:9

   indx = find(dat_sum==i);
   forhist(i) = length(indx);

end forhist = forhist';

% ITK Snap Label colors mycolors = [255 2 0; 228 96 0; 231 178 0; 224 255 10; 87 255 21; 0 255 145; 35 252 255; 7 81 255; 85 36 255]./255;

% Histogram figure(1), h = bar(forhist) colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('number of NAc voxels') title('Subject overlap in NAc voxels (1 x 1 x 1 mm)') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(forhist,1); for j = 1:length(forhist)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)

clear h ch fvd fvcd zs izs j

% Cumulative Histogram c_elements = cumsum(forhist); figure(2), h = bar(c_elements) colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('number of NAc voxels') title('Cumulative Histogram of subject overlap in NAc voxels (1 x 1 x 1 mm)') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(c_elements,1); for j = 1:length(c_elements)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)

clear h ch fvd fvcd zs izs j

% proportional histogram (rhist) rhist_bars = forhist./sum(forhist); figure(3), h = bar(rhist_bars); colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('proportion of total NAc voxels') title('Relative histogram of subject overlap in NAc voxels (1 x 1 x 1 mm)') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(rhist_bars,1); for j = 1:length(rhist_bars)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)

clear h ch fvd fvcd zs izs j c_elements forhist i indx indx2 k temp

%% 2x2x2 voxels

% Size - finds the # of voxels in each subjects mask then calculates the average and the s.e.

for k = 1:size(dat2,4);  % 2x2x2

   indx = find(dat2(:,:,:,k));
   subj_vox_count2(k) = length(indx);
   

end

mean_nvox2 = mean(subj_vox_count2); se_nvox2 = std(subj_vox_count2)/(size(dat2,4)^.5);


%% PLOTS

forhist = [];

for i=1:9

   indx = find(dat_sum==i);
   forhist(i) = length(indx);

end forhist = forhist';

% ITK Snap Label colors mycolors = [255 2 0; 228 96 0; 231 178 0; 224 255 10; 87 255 21; 0 255 145; 35 252 255; 7 81 255; 85 36 255]./255;

% Histogram figure(4), h = bar(forhist) colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('number of NAc voxels') title('Subject overlap in NAc voxels (2 x 2 x 2 mm)') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(forhist,1); for j = 1:length(forhist)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)

clear h ch fvd fvcd zs izs j

% Cumulative Histogram c_elements = cumsum(forhist); figure(5), h = bar(c_elements) colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('number of NAc voxels') title('Cumulative Histogram of subject overlap in NAc voxels (2 x 2 x 2 mm)') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(c_elements,1); for j = 1:length(c_elements)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)

clear h ch fvd fvcd zs izs j

% proportional histogram (rhist) rhist_bars2 = forhist./sum(forhist); figure(6), h = bar(rhist_bars2); colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('proportion of total NAc voxels') title('Relative histogram of subject overlap in NAc voxels (2 x 2 x 2 mm)') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(rhist_bars2,1); for j = 1:length(rhist_bars2)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)


%% Compare to single subject rep NAc

in_indx = find(dat_rep);  % # of NAc voxels in Colin mask in_rep_mask = dat_sum(in_indx); out_indx = find(dat_rep==0); out_temp = dat_sum(out_indx); out_indx2 = find(out_temp);  % # of NAc voxels out of Colin mask out_rep_mask = out_temp(out_indx2);

mean_in_mask = mean(in_rep_mask);  % mean subject overlap/NAc voxel within Colin mask mean_out_mask = mean(out_rep_mask); % mean subject overlap/NAc voxel outside of mask


clear h ch fvd fvcd zs izs j c_elements forhist i indx indx2 k temp

% plot relative hist for reg and double size voxels next to each other

% figure(7) % ax(1) = subplot(1,2,1); % rgb = imread('ngc6543a.jpg'); % image(rgb); title('RGB image') % ax(2) = subplot(122); % im = mean(rgb,3); % image(im); title('Intensity Heat Map') % colormap(hot(256)) % linkaxes(ax,'xy') % axis(ax,'image')

% proportional histogram (rhist)

figure(7) subplot(121), h = bar(rhist_bars); colormap(mycolors); xlabel('number of subjects overlap/ NAc voxel') ylabel('proportion of total NAc voxels') title('1 x 1 x 1 mm voxels ') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(rhist_bars,1); for j = 1:length(rhist_bars)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd) subplot(122),h = bar(rhist_bars2); colormap(mycolors); title('2 x 2 x 2 mm voxels') ch = get(h,'Children'); %getting a handle on the kids fvd = get(ch,'Faces'); fvcd = get(ch,'FaceVertexCData'); [zs, izs] = sortrows(rhist_bars2,1); for j = 1:length(rhist_bars2)

   row = izs(j);
   fvcd(fvd(row,:)) = j;

end set(ch,'FaceVertexCData',fvcd)