bfm_plot_datamatcorrs_3v.m

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

   brainlv_lst = {};

   for g = 1:num_grp

      brainlv = datamatcorrs_lst{g};		% borrow 'blv' name to reuse code

      r = size(brainlv, 1) / num_conditions;
      c = size(brainlv, 2);

      for b = 1:num_behav

         mask = [0:(num_conditions-1)]*num_behav + b;
         tmp = brainlv(mask, :);
         tmp = tmp';
         brainlv_lst{g,b} = tmp;

%      for i = 1:num_conditions
%         tmp{i} = brainlv(r*(i-1)+1:r*i, :);
%      end  
%
%      brainlv = ones(num_conditions, c);
%
%      for i = 1:num_conditions
%         brainlv(i,:) = mean(tmp{i},1);
%      end
%
%      brainlv = brainlv';
%      brainlv_lst{g} = brainlv;

      end
   end


   setappdata(gcf,'brainlv',brainlv_lst);
   setappdata(gcf, 's', s);
   set_blv_fields(grp_idx,behav_idx,lv_idx);

   if ~exist('boot_result','var') | isempty(boot_result)
      ToggleView(0);
      set(findobj(gcf,'Tag','ViewMenu'),'Visible','off');
   else					% show bootstrap ratio if exist
      ToggleView(1);
      set(findobj(gcf,'Tag','ViewMenu'),'Visible','on');

      % set the bootstrap ratio field values
      %
      setappdata(gcf,'BSRatio',boot_result.compare);
      set_bs_fields(bs_lv_idx);
      UpdatePValue;
   end;

   h = findobj(gcf,'Tag','OpenContrastWindow');
   if 1 % isequal(ContrastFile,'NONE') | isequal(ContrastFile,'BEHAV')
      set(h,'Visible','off');
   else
      set(h,'Visible','on');
   end;

   set(gcf,'Pointer',old_pointer);
   set(findobj(gcf,'Tag','MessageLine'),'String','');

   setappdata(gcf,'SessionFileList', SessionProfiles);
   setappdata(gcf,'RotateAmount',rot_amount);
   setappdata(gcf,'CurrGroupIdx',grp_idx);
   setappdata(gcf,'CurrBehavIdx',behav_idx);
   setappdata(gcf,'NumGroup',num_grp);
   setappdata(gcf,'CurrLVIdx',lv_idx);
   setappdata(gcf,'CurrBSLVIdx',bs_lv_idx);
   setappdata(gcf,'STDims',st_dims);
   setappdata(gcf,'NumBSLVs',bs_num_lv)

   return;						% load_pls_result


%-------------------------------------------------------------------------
%
function OpenResponseFnPlot()

  sessionFileList = getappdata(gcbf,'SessionFileList');

  rf_plot = getappdata(gcf,'RFPlotHdl');
  if ~isempty(rf_plot)
      msg = 'ERROR: Response function plot is already been opened';
      set(findobj(gcf,'Tag','MessageLine'),'String',msg);
      return;
  end;

  rf_plot = fmri_plot_rf('LINK',sessionFileList);
  link_info.hdl = gcbf;
  link_info.name = 'RFPlotHdl';
  setappdata(rf_plot,'LinkFigureInfo',link_info);
  setappdata(gcbf,'RFPlotHdl',rf_plot);

  %  make sure the Coord of the Response Function Plot contains 
  %  the current point in the Response
  %
  cur_coord = getappdata(gcf,'Coord');
  setappdata(rf_plot,'Coord',cur_coord);

  return;					% OpenResponseFnPlot


%-------------------------------------------------------------------------
%
function SetClusterReportOptions()

   st_origin = getappdata(gcbf,'Origin');
   st_dims = getappdata(gcbf,'STDims');
   cluster_min_size = getappdata(gcbf,'ClusterMinSize');
   cluster_min_dist = getappdata(gcbf,'ClusterMinDist');

   if isempty(st_origin) | all(st_origin == 0)

      st_voxel_size = getappdata(gcf,'STVoxelSize');

      if all(st_dims == [40 48 1 34]) & all(st_voxel_size == [4 4 4])
         st_origin = [20 29 12];
      elseif all(st_dims == [91 109 1 91]) & all(st_voxel_size == [2 2 2])
         st_origin = [46 64 37];
      else
         % according to SPM: if the origin field contains 0, then the origin is
         % assumed to be at the center of the volume.
         %
         st_origin = floor((dims([1 2 4])+1)/2);
         % st_origin = round(st_dims([1 2 4])/2);
      end;
   end;

   prompt = {'Minimum cluster size (in voxels)',  ...
	     'Minimum distance (in mm) between cluster peaks', ...
	     'Origin location (in voxels)' };
   defValues = { num2str(cluster_min_size), ...
		 num2str(cluster_min_dist), ...
		 num2str(st_origin)};
   dlgTitle='Cluster Report Options';
   lineNo = 1;
   answer = inputdlg(prompt,dlgTitle,lineNo,defValues);

   if isempty(answer),
      return;
   end;

   invalid_options = 0;
   min_size = str2num(answer{1}); 
   min_dist = str2num(answer{2}); 
   origin_xyz = str2num(answer{3}); 

   if isempty(min_size) | isempty(min_dist) | isempty(origin_xyz)
      invalid_options = 1;
   elseif (min_size <= 0) | (min_dist <= 0) | (sum(origin_xyz<= 0) ~= 0)
      invalid_options = 1;
   end;   
   
   if (invalid_options)
	msg = 'Invalid cluster report options.  Options do not changed';
	set(findobj(gcbf,'Tag','MessageLine'),'String',msg);
	return;
   end;

   setappdata(gcbf,'Origin',origin_xyz);
   setappdata(gcbf,'ClusterMinSize',min_size);
   setappdata(gcbf,'ClusterMinDist',min_dist);

   return;					% SetClusterReportOptions


%-------------------------------------------------------------------------
%
function OpenClusterReport()

   %  wait message
   old_pointer = get(gcbf,'Pointer');
   set(gcbf,'Pointer','watch');

   msg = 'Generating Cluster Report ... please wait';
   h = rri_wait_box(msg, [0.5 0.1]);

   cluster_min_size = getappdata(gcbf,'ClusterMinSize');
   cluster_min_dist = getappdata(gcbf,'ClusterMinDist');

   fmri_cluster_report(cluster_min_size,cluster_min_dist);


   set(gcbf,'Pointer',old_pointer);
   set(findobj(gcbf,'Tag','MessageLine'),'String','');

   delete(h);

   return;					% OpenClusterReport


%-------------------------------------------------------------------------
function OpenContrastWindow()

   contrast_fig = getappdata(gcbf,'ContrastFigHdl');
   if ~isempty(contrast_fig)
      msg = 'ERROR: Constrasts information has already been dispalyed.';
      set(findobj(gcf,'Tag','MessageLine'),'String',msg);
      return;
   end  

   h = findobj(gcbf,'Tag','ResultFile');
   PLSresultFile = get(h,'UserData');

   load(PLSresultFile,'ContrastFile','SessionProfiles');

   if isequal(ContrastFile,'NONE'), 
      msg = 'No contrast was used for this PLS analysis.'; 
      set(findobj(gcf,'Tag','MessageLine'),'String',msg);
      return;
   end;

   if isequal(ContrastFile,'HELMERT'),   % using Helmert matrix for contrasts
      load(SessionProfiles{1}{1});

      conditions = session_info.condition;
      num_conditions = length(conditions);
      helmert_contrasts = rri_helmert_matrix(num_conditions);

      for i=1:num_conditions-1,
         pls_contrasts(i).name = sprintf('Contrast #%d',i);
         pls_contrasts(i).value = helmert_contrasts(:,i)';
      end;
   else
      try
         load(ContrastFile);
      catch 
         msg = sprintf('ERROR: Cannot open contrast file "%s".',ContrastFile); 
         set(findobj(gcf,'Tag','MessageLine'),'String',msg);
         return;
      end;
   end;

   contrast_fig = fmri_input_contrast_ui(pls_contrasts,conditions,1);

   link_info.hdl = gcbf;
   link_info.name = 'ContrastFigHdl';
   setappdata(contrast_fig,'LinkFigureInfo',link_info);
   setappdata(gcbf,'ContrastFigHdl',contrast_fig);

   return;					% OpenContrastWindow


%-------------------------------------------------------------------------
%
function OpenDesignPlot()


   scores_fig = getappdata(gcbf,'ScorePlotHdl');
   if ~isempty(scores_fig)
      msg = 'ERROR: Design Scores Plot is already been opened';
      set(findobj(gcf,'Tag','MessageLine'),'String',msg);
      return;
   end  

   h = findobj(gcbf,'Tag','ResultFile');
   PLSresultFile = get(h,'UserData');

   scores_fig = fmri_plot_scores('LINK',PLSresultFile);

   lv_idx = getappdata(gcbf,'CurrLVIdx');
   if (lv_idx ~= 1)
      fmri_plot_scores('UPDATE_LV_SELECTION',scores_fig,lv_idx);
   end;

   link_info.hdl = gcbf;
   link_info.name = 'ScorePlotHdl';
   setappdata(scores_fig,'LinkFigureInfo',link_info);
   setappdata(gcbf,'ScorePlotHdl',scores_fig);

   return;					% OpenDesignPlot


%-------------------------------------------------------------------------
%
function OpenBrainScoresPlot()


  bs_plot = getappdata(gcf,'BSPlotHdl');
  if ~isempty(bs_plot)
      msg = 'ERROR: Brain score plot is already been opened';
      set(findobj(gcf,'Tag','MessageLine'),'String',msg);
      return;
  end;

  sessionFileList = getappdata(gcbf,'SessionFileList');

  h = findobj(gcbf,'Tag','ResultFile');
  PLSresultFile = get(h,'UserData');


  bs_plot = fmri_plot_brain_scores('LINK',sessionFileList,PLSresultFile);
  link_info.hdl = gcbf;
  link_info.name = 'BSPlotHdl';
  setappdata(bs_plot,'LinkFigureInfo',link_info);
  setappdata(gcbf,'BSPlotHdl',bs_plot);

  return;					% OpenBrainScoresPlot


%-------------------------------------------------------------------------
%
function SaveResultToIMG(is_disp)
%

  h = findobj(gcf,'Tag','ResultFile'); PLSresultFile = get(h,'Userdata');

  try 				% load the dimension info of the st_datamat
     load(PLSresultFile,'st_dims'),
  catch
     msg =sprintf('ERROR: Cannot load the PLS result file "%s".',PLSresultFile);
     set(findobj(gcf,'Tag','MessageLine'),'String',msg);
     return;
  end;

  h = findobj(gcf,'Tag','LVIndexEdit');  lv_idx = get(h,'Userdata');
  curr_lv_idx = getappdata(gcf,'CurrLVIdx');
  if (lv_idx ~= curr_lv_idx),
     lv_idx = curr_lv_idx;
     set(h,'String',num2str(lv_idx));
  end;

  h = findobj(gcf,'Tag','GroupIndexEdit');  grp_idx = get(h,'Userdata');
  curr_grp_idx = getappdata(gcf,'CurrGroupIdx');
  if (grp_idx ~= curr_grp_idx),
     grp_idx = curr_grp_idx;
     set(h,'String',num2str(grp_idx));
  end;

  h = findobj(gcf,'Tag','BehavIndexEdit');  behav_idx = get(h,'Userdata');
  curr_behav_idx = getappdata(gcf,'CurrBehavIdx');
  if (behav_idx ~= curr_behav_idx),
     behav_idx = curr_behav_idx;
     set(h,'String',num2str(behav_idx));
  end;

  old_pointer = get(gcf,'Pointer');
  fig_hdl = gcf;
  set(fig_hdl,'Pointer','watch');

  if 1 % (getappdata(gcf,'ViewBootstrapRatio') == 0),	% save brain lv 
     thresh = getappdata(gcf,'BLVThreshold'); 
     if is_disp
        create_st_brainlv_disp(PLSresultFile,lv_idx,thresh,grp_idx,behav_idx);
     else
        create_st_brainlv_img(PLSresultFile,lv_idx,thresh,grp_idx,behav_idx);
     end
  else							% save bootstrap ratio
     thresh_ratio = getappdata(gcf,'BSThreshold');  
     create_bs_ratio_img(PLSresultFile,lv_idx,thresh_ratio);
  end;

  set(fig_hdl,'Pointer',old_pointer);

  return;					% SaveResultToIMG


%--------------------------------------------------------------------------
function create_st_brainlv_disp(PLSresultFile,lv_idx,thresh_ratio,grp_idx,behav_idx);

  %  get the output IMG filename first  
  %
  [pn fn] = fileparts(PLSresultFile);
  resultfile_prefix = fn(1:end-11);

  image_fn = sprintf('%sBfMRIdatcorr_disp_grp%d_cond%d_beh%d.img',resultfile_prefix,grp_idx,lv_idx,behav_idx);

  [filename, pathname] = uiputfile(image_fn,'Datamat Correlation IMG file');

  img_file = [pathname, filesep, filename];

  if isequal(filename,0)
      return;
  end;

  %  load the result file
  %
  load(PLSresultFile,'st_dims','st_win_size','st_coords', ...
		     'st_voxel_size','st_origin');

  dims = st_dims([1 2 4]);
  newcoords = st_coords;

  blv = getappdata(gcbf,'BLVData');
  brainlv = blv{grp_idx, behav_idx};

  %  save background to img
  %
  brainlv = brainlv(:,lv_idx);


   bs = getappdata(gcbf,'BSRatio');
   h = findobj(gcf,'Tag','BSLVIndexEdit'); bs_lv_idx = str2num(get(h,'String'));
   h = findobj(gcf,'Tag','BSThreshold'); bs_thresh = str2num(get(h,'String'));
   bs = bs(:, bs_lv_idx);
   bs_strong = zeros(size(bs));
   bs_idx = [find(bs <=- bs_thresh); find(bs >= bs_thresh)];
   bs_strong(bs_idx) = 1;
   brainlv = brainlv .* bs_strong;


  bg_img = getappdata(gcf,'BackgroundImg');
  cmap = getappdata(gcf,'cmap');
  num_blv_colors = 25;
  brain_region_color_idx = 51;
  first_lower_color_idx = 101;
  first_upper_color_idx = 126;
  h = findobj(gcf,'Tag','MaxValue'); max_blv = str2num(get(h,'String'));
  h = findobj(gcf,'Tag','MinValue'); min_blv = str2num(get(h,'String'));
  too_large = find(brainlv > max_blv); brainlv(too_large) = max_blv;
  too_small = find(brainlv < min_blv); brainlv(too_small) = min_blv;

  % Create the image slices in which voxels are set to be within certain range
  %
  lower_interval = (abs(min_blv) - thresh_ratio) / (num_blv_colors-1);
  upper_interval = (max_blv - thresh_ratio) / (num_blv_colors-1);

  blv = zeros(1,length(newcoords)) + brain_region_color_idx;
  lower_idx = find(brainlv <= -thresh_ratio);
  blv_offset = brainlv(lower_idx) - min_blv; 
  lower_color_idx = round(blv_offset/lower_interval)+first_lower_color_idx;
  blv(lower_idx) = lower_color_idx;

  upper_idx = find(brainlv >= thresh_ratio);
  blv_offset = max_blv - brainlv(upper_idx); 
  upper_color_idx = num_blv_colors - round(blv_offset/upper_interval);
  upper_color_idx = upper_color_idx + first_upper_color_idx - 1;
  blv(upper_idx) = upper_color_idx;

  if isempty(bg_img)
     non_brain_region_color_idx = size(cmap,1);
     img = zeros(1,dims(1)*dims(2)*dims(3)) + non_brain_region_color_idx;
     img(newcoords) = blv;
     img = reshape(img,dims); 
  else
     max_bg = max(bg_img(:));
     min_bg = min(bg_img(:));
     img = (bg_img - min_bg) / (max_bg - min_bg) * 100;
     img(newcoords(lower_idx)) = blv(lower_idx);
     img(newcoords(upper_idx)) = blv(upper_idx);
     img = reshape(img,dims); 
  end

%  img = zeros(dims);

%  blv = brainlv(1,lv_idx); 
%  blv(abs(blv) < thresh_ratio) = 0;

%  img(st_coords) = blv;
     
  descrip = sprintf('DatamatCorrelation from %s, Group: %d, Condition: %d, Behavior: %d, Threshold: %8.5f', ...
		PLSresultFile,grp_idx,lv_idx,behav_idx,thresh_ratio);
  rri_write_img(img_file,img,0,dims,st_voxel_size,16,st_origin,descrip);

  %  save background to img
  %
  [tmp filename] = fileparts(filename);
  save(filename,'cmap');

  return;					% create_st_brainlv_disp