bfm_plot_rf_task.m

绝对经典,老外制作的功能强大的matlab实现PLS_TOOBOX

      else
         cluster_info = getappdata(CallBackFig, 'cluster_blv');
      end

      %  Get cluster voxels coords relative to whole volume
      %
      if length(cluster_info) < CurrLVIdx
         cluster_info = [];
      else
         cluster_info = cluster_info{CurrLVIdx};
         cluster_info = cluster_info.data{1}.idx;
      end

      %  If "Bootstrap" is computed, voxels that meet the bootstrap ratio
      %  threshold will be used as a criteria to select surrounding voxels
      %
      BSThreshold = getappdata(CallBackFig,'BSThreshold');

      if ~isempty(BSThreshold)
         BSRatio = getappdata(CallBackFig,'BSRatio');
         BSRatio = BSRatio(:,CurrLVIdx);

         all_voxel = zeros(1,prod(dims));
         all_voxel(common_coords) = (abs(BSRatio) >= BSThreshold);
         bsridx = find(all_voxel);
      else
         bsridx = [];
      end

      %  Only including surrounding voxels that meet the bootstrap ratio
      %  threshold, or are part of cluster masked voxels
      %
      bsr_cluster_coords = unique([cluster_info bsridx cur_coord]);

      %  Intersect of neighborhood coord "neighbor_coord" and 
      %  "bsr_cluster_coords"
      %
      neighbor_coord = intersect(neighbor_coord, bsr_cluster_coords);

   else

      neighbor_coord = cur_coord;

   end;		% if neighbor_size > 0

   %  find out neighborhood indices in st_datamat
   %
   [ncoord ncoord_idx] = intersect(common_coords, neighbor_coord);

   neighbor_numbers = length(ncoord_idx);
   h = findobj(gcf,'Tag','neighbornumberEdit');
   set(h,'String',num2str(neighbor_numbers));

   st_datamat = getappdata(gcf,'STDatamat'); 
   st_data = st_datamat;
   selected_condition = getappdata(gcf,'SelectedCondition');

   % generate the plots
   %
   cond_idx = find(selected_condition == 1);

   setappdata(gcf,'ST_data',st_data);
   setappdata(gcf,'PlotCondIdx',cond_idx);

   if isempty(getappdata(gcf,'CombinePlots'))
      setappdata(gcf,'CombinePlots',0);
   end;

   if isempty(getappdata(gcf,'ShowAverage'))
      setappdata(gcf,'ShowAverage',0);
   end;

   %  the following code is to get an intensity array for the voxel and plot it
   %
   num_subj = getappdata(gcf,'num_subj');
   subj_name = getappdata(gcf,'subj_name');
   subjects = [subj_name, {'mean'}];

   for k = cond_idx			% cond
      for n = 1:num_subj		% subj
         j = n+(k-1)*num_subj;		% row number in datamat
         intensity(k,n) = mean(st_data(j, ncoord_idx), 2);
         intensity_hdl(k,n) = 0;	% initialization
      end
      intensity_avg_hdl(k) = 0;
   end

   intensity_avg = mean(intensity,2);

   color_code =[ 'bo';'rx';'g+';'m*';'bs';'rd';'g^';'m<';'bp';'r>'; ...
                 'gh';'mv';'ro';'gx';'m+';'b*';'rs';'gd';'m^';'b<'];

   % need more color
   %
   if num_subj+1 > size(color_code,1)

      tmp = [];

      for i=1:ceil((num_subj+1)/size(color_code,1))
         tmp = [tmp; color_code];
      end

      color_code = tmp;

   end

   cla; grid off; hold on;
   for k = cond_idx

      for n = 1:num_subj
         intensity_hdl(k,n) = plot(k,intensity(k,n), ...
			color_code(n,:));
      end

      intensity_avg_hdl(k) = bar(k,intensity_avg(k));
      set(intensity_avg_hdl(k),'facecolor','none')
      % set(intensity_avg_hdl(k), 'linewidth', 2);

   end

   % normalize with intensity(:), not st_data(:)
   %
   min_y = min(intensity(:)); max_y = max(intensity(:));
   margin_y = abs((max_y - min_y) / 100);

   axis_scale = getappdata(gcf, 'axis_scale');

   if isempty(axis_scale)
      axis_scale = [0 length(cond_idx)+1 min_y-margin_y max_y+margin_y];
   else
      if (min_y-margin_y)<axis_scale(3)
         axis_scale(3) = min_y-margin_y;
      end;

      if (max_y+margin_y)>axis_scale(4)
         axis_scale(4) = max_y+margin_y;
      end;
   end

   setappdata(gcf, 'axis_scale', axis_scale);
   axis(axis_scale);

   set(gca, 'xtick', cond_idx);

   xlabel('Conditions');
   ylabel('Intensities');

   xyz = getappdata(gcf,'XYZ');
   title(['Intensity for subjects and mean at voxel:  [',num2str(xyz),']']);

   hold off;

   l_hdl = [];

   if ~isempty(subjects),

      intensity_legend = [intensity_hdl(1,:), intensity_avg_hdl(1)];

      % remove the old legend to avoid the bug in the MATLAB5
      old_legend = getappdata(gcf,'LegendHdl');
      if ~isempty(old_legend),
        try
          delete(old_legend{1});
        catch
        end;
      end;

      % create a new legend, and save the handles
      [l_hdl, o_hdl] = legend(intensity_legend, subjects, 0);
      legend_txt(o_hdl);
      set(l_hdl,'color',[0.9 1 0.9]);
      setappdata(gcf,'LegendHdl',[{l_hdl} {o_hdl}]);

   else

      setappdata(gcf,'LegendHdl',[]);

   end;

   setappdata(gcf, 'ncoord_idx', ncoord_idx);
   setappdata(gcf, 'ncoord', ncoord);
   set(findobj(gcf,'Tag','DataMenu'),'visible','on');

   return;					% plot_response_fn


%--------------------------------------------------------------------------
%
function  save_response_fn()

   neighbor_numbers = 1;

   %  Get Neighborhood Size
   h = findobj(gcf,'Tag','neighborhoodEdit');
   neighbor_size = round(str2num(get(h,'String')));

   if isempty(neighbor_size) | ~isnumeric(neighbor_size)
      neighbor_size = 0;
   end

   st_coords = getappdata(gcf,'STCoords');
   common_coords = getappdata(gcf,'common_coords');

   xyz = getappdata(gcf,'XYZ');

   %  extract the currect ploting data
   %
   cur_coord = getappdata(gcf,'Coord');
   if isempty(cur_coord)
       msg = 'ERROR: No point has been seleted to plot.';
       set(findobj(gcf,'Tag','MessageLine'),'String',msg);
       return;
   end;

   coord_idx = find(common_coords == cur_coord);
   if isempty(coord_idx)
       msg = 'ERROR: The selected point is outside the brain region.';
       set(findobj(gcf,'Tag','MessageLine'),'String',msg);
       return;
   end;
%   coord_idx = find(st_coords == cur_coord);
   ncoord_idx = getappdata(gcf, 'ncoord_idx');
   neighbor_numbers = length(ncoord_idx);

   CallBackFig = getappdata(gcf,'CallBackFig');
   dims = getappdata(CallBackFig,'STDims');
   ncoord = getappdata(gcf, 'ncoord');
   xyz = rri_coord2xyz(ncoord,dims);

   %  get selected file names 
   %
   st_files = getappdata(gcf,'DatamatFile');
   h = findobj(gcf,'Tag','STDatamatPopup');
   selected_idx = get(h,'Value');
   selected_files = st_files{selected_idx};

   num_cond = getappdata(gcf,'num_cond');
   num_subj_lst = getappdata(gcf,'num_subj_lst');
   num_subj = num_subj_lst(selected_idx);

   mask = 1:num_subj*num_cond;
   mask = reshape(mask, [num_subj, num_cond]);

   %  extract the data of the selected voxel
   %
   st_data = getappdata(gcf,'STDatamat');
   data = mean(st_data(:, ncoord_idx), 2);
   standard_deviation = std(st_data(:, ncoord_idx), 0, 2);

   %  save behavdata
   %
   for i = 1:size(mask,1)
      behavdata.subject{i} = [data(mask(i,:))]';
   end;

   [filename, pathname] = ...
	uiputfile('*_BfMRI_rf_plot.mat','Save the Response Functions');

   if ischar(filename) & isempty(findstr(lower(filename),'_bfmri_rf_plot'))
      [tmp filename] = fileparts(filename);
      filename = [filename, '_BfMRI_rf_plot.mat'];
   end

   if isequal(filename,0)
	return;
   end;
   rf_plot_file = fullfile(pathname,filename);

   xyz_str = 'xyz';
 
   sa = getappdata(gcf, 'sa');
   
   if ~isempty(sa) & sa == 1
      yzx = xyz;
      xyz_str = 'yzx';
   elseif ~isempty(sa) & sa == 0
      xzy = xyz;
      xyz_str = 'xzy';
   end

   try
     save (rf_plot_file, 'selected_files', 'data', 'standard_deviation', 'xyz', 'neighbor_size', 'neighbor_numbers' );
   catch
     msg = sprintf('Cannot save the response function data to %s',rf_plot_file);
     set(findobj(gcf,'Tag','MessageLine'),'String',['ERROR: ', msg]);
     status = 0;
     return;
   end;

   return;                                      % save_response_fn


%--------------------------------------------------------------------------
%
function  save_response_fn2()

   st_coords = getappdata(gcf,'STCoords');
   common_coords = getappdata(gcf,'common_coords');

   xyz = getappdata(gcf,'XYZ');

   %  extract the currect ploting data
   %
   cur_coord = getappdata(gcf,'Coord');
   if isempty(cur_coord)
       msg = 'ERROR: No point has been seleted to plot.';
       set(findobj(gcf,'Tag','MessageLine'),'String',msg);
       return;
   end;

   coord_idx = find(common_coords == cur_coord);
   if isempty(coord_idx)
       msg = 'ERROR: The selected point is outside the brain region.';
       set(findobj(gcf,'Tag','MessageLine'),'String',msg);
       return;
   end;
%   coord_idx = find(st_coords == cur_coord);

   %  get selected file names 
   %
   st_files = getappdata(gcf,'DatamatFile');
   h = findobj(gcf,'Tag','STDatamatPopup');
   selected_idx = get(h,'Value');
   selected_files = st_files{selected_idx};

   num_cond = getappdata(gcf,'num_cond');
   num_subj_lst = getappdata(gcf,'num_subj_lst');
   num_subj = num_subj_lst(selected_idx);

   mask = 1:num_subj*num_cond;
   mask = reshape(mask, [num_subj, num_cond]);

   %  extract the data of the selected voxel
   %
   st_data = getappdata(gcf,'STDatamat');
   data = st_data(:,coord_idx);

   %  save behavdata
   %
   for i = 1:size(mask,1)
      behavdata.subject{i} = [data(mask(i,:))]';
   end;

   pattern = ...
      ['<ONLY INPUT PREFIX>*_BfMRI_grp',num2str(selected_idx),'_subj1_behavdata1.txt'];

   [filename, pathname] = uiputfile(pattern,'Save Behav Data');

   if isequal(filename,0)
	return;
   end;

   [tmp filename] = fileparts(filename);

   for i = 1:size(mask,1)
      rf_plot_file = fullfile(pathname, ...
	sprintf('%s_BfMRI_grp%d_subj%d_behavdata1.txt',filename,selected_idx,i));
      behavdata = double(data(mask(i,:)));
      save (rf_plot_file, '-ascii', 'behavdata');
   end

   return;                                      % save_response_fn2


%--------------------------------------------------------------------------
%
function  save_response_behav()

   xyz = getappdata(gcf,'XYZ');

   %  extract the currect ploting data
   %
   cur_coord = getappdata(gcf,'Coord');
   if isempty(cur_coord)
       msg = 'ERROR: No point has been seleted to plot.';
       set(findobj(gcf,'Tag','MessageLine'),'String',msg);
       return;
   end;

   %  find out the coord_idx
   %
   st_coords = getappdata(gcf,'STCoords');
   common_coords = getappdata(gcf,'common_coords');

   coord_idx = find(common_coords == cur_coord);
   if isempty(coord_idx)
       msg = 'ERROR: The selected point is outside the brain region.';
       set(findobj(gcf,'Tag','MessageLine'),'String',msg);
       return;
   end;
   coord_idx = find(st_coords == cur_coord);

   st_data = getappdata(gcf,'datamatcorrs');
   data = st_data(:,coord_idx);

   [filename, pathname] = ...
	uiputfile('*_BfMRI_datamatcorr_plot.mat','Save the Response Functions');

   if ischar(filename) & isempty(findstr(lower(filename),'_bfmri_datamatcorr_plot'))
      [tmp filename] = fileparts(filename);
      filename = [filename, '_BfMRI_datamatcorr_plot.mat'];
   end

   if isequal(filename,0)
	return;
   end;
   rf_plot_file = fullfile(pathname,filename);

   try
     save (rf_plot_file, 'data', 'xyz' );
   catch
     msg = sprintf('Cannot save the response function data to %s',rf_plot_file);
     set(findobj(gcf,'Tag','MessageLine'),'String',['ERROR: ', msg]);
     status = 0;
     return;
   end;

   return;                                      % save_response_behav


%--------------------------------------------------------------------------
%
function   set_cond_axes(num_rows,num_cols,axes_margin)
%
%   Define the axes for the response functions of different conditions
%
%   axes_margin: [left right bottom top], which specified in terms of 
%                  normal unit respect to the figure;
% 

  show_avg = getappdata(gcf,'ShowAverage');
  if ~isempty(show_avg) & (show_avg == 1),
     num_cols = num_cols + 1;
  end;

  %  clear up the old handles
  %
  old_hdls = getappdata(gcf,'AxesHdls');
  if ~isempty(old_hdls)
     for i=1:length(old_hdls(:))
	 if ishandle(old_hdls{i}), delete(old_hdls{i}); end
     end; 
  end;

  old_hdls = getappdata(gcf,'AxesCombineHdls');
  if ~isempty(old_hdls)
     for i=1:length(old_hdls(:))
	 if ishandle(old_hdls{i}), delete(old_hdls{i}); end
     end; 
  end;

  f_pos = get(gcf,'Position');
  axes_boundary(1) = axes_margin(1);
  axes_boundary(2) = 1 - axes_margin(2);
  axes_boundary(3) = axes_margin(3);
  axes_boundary(4) = 1 - axes_margin(4);

  %  plot data in each axis
  %
  ax_hdls = cell(num_rows,num_cols);
  ax_pos = cell(num_rows,num_cols);
  ax_combine_hdls = cell(num_rows,1);
  ax_combine_pos = cell(num_rows,1);

  plot_width  = (axes_boundary(2) - axes_boundary(1)) / num_cols;
  plot_height = (axes_boundary(4) - axes_boundary(3)) / num_rows;
  axes_min_x = axes_boundary(1);
  axes_min_y = axes_boundary(3);

  for row=1:num_rows, 

    axes_y = axes_min_y + plot_height*(num_rows-row); 

    % for separate plots within each condition 
    %
    for col=1:num_cols,