CellsortChoosePCs

 [PCuse] = CellsortChoosePCs(fn, mixedfilters)

 Allows the user to select which principal components will be kept
 following dimensional reduction.

 Inputs:
   fn - movie file name. Must be in TIFF format.
   mixedfilters - N x X matrix of N spatial signal mixtures sampled at X
   spatial points.

 Outputs:
   PCuse - vector of indices of the PCs to be kept for dimensional
   reduction

 Eran Mukamel, Axel Nimmerjahn and Mark Schnitzer, 2009
 Email: eran@post.harvard.edu, mschnitz@stanford.edu

0001 function [PCuse] = CellsortChoosePCs(fn, mixedfilters)
0002 % [PCuse] = CellsortChoosePCs(fn, mixedfilters)
0003 %
0004 % Allows the user to select which principal components will be kept
0005 % following dimensional reduction.
0006 %
0007 % Inputs:
0008 %   fn - movie file name. Must be in TIFF format.
0009 %   mixedfilters - N x X matrix of N spatial signal mixtures sampled at X
0010 %   spatial points.
0011 %
0012 % Outputs:
0013 %   PCuse - vector of indices of the PCs to be kept for dimensional
0014 %   reduction
0015 %
0016 % Eran Mukamel, Axel Nimmerjahn and Mark Schnitzer, 2009
0017 % Email: eran@post.harvard.edu, mschnitz@stanford.edu
0018 
0019 fprintf('-------------- CellsortICA %s -------------- \n', date)
0020 
0021 [pixw,pixh] = size(imread(fn,1));
0022 
0023 npcs = 20; % Number of PCs to display concurrently
0024 currpcs = [1:npcs];
0025 PCf = [];
0026 while isempty(PCf)
0027     showpcs(currpcs, mixedfilters, pixw, pixh)
0028     yl = ylim;
0029     xl = xlim;
0030     set(gca,'Units','pixels')
0031     title(['Choose first PC; showing PCs [',num2str(currpcs(1)),':',num2str(currpcs(end)),']'])
0032     PCf = input('Number of first PC to retain, Klow (''b/f'' to scroll backwards/forwards)): ','s');
0033     if PCf=='b'
0034         currpcs = currpcs - min(npcs,currpcs(1)-1);
0035         PCf = [];
0036     elseif (PCf=='f')
0037         currpcs = currpcs+npcs;
0038         if nnz(currpcs>size(mixedfilters,2))
0039             currpcs = [-npcs+1:0]+size(mixedfilters,2);
0040             fprintf('Reached end of stored PCs.\n')
0041         end
0042         PCf = [];
0043     else
0044         PCf = str2num(PCf);
0045     end
0046 end
0047 PCl=[];
0048 currpcs = [PCf:PCf+npcs-1];
0049 while isempty(PCl)
0050     showpcs(currpcs, mixedfilters, pixw, pixh)
0051     title(['Choose last PC; showing PCs [',num2str(currpcs(1)),':',num2str(currpcs(end)),']'])
0052     PCl = input('Number of last PC to retain, Khigh (''b/f'' to scroll backwards/forwards): ','s');
0053     if PCl=='b'
0054         currpcs = currpcs - min(npcs,currpcs(1)-1);
0055         PCl = [];
0056     elseif (PCl=='f')
0057         currpcs = currpcs+npcs;
0058         if nnz(currpcs>size(mixedfilters,2))
0059             currpcs = [-npcs+1:0]+size(mixedfilters,2);
0060             fprintf('Reached end of stored PCs.\n')
0061         end
0062         PCl = [];
0063     else
0064         PCl = str2num(PCl);
0065     end
0066 end
0067 currpcs = [PCf:PCl];
0068 PCbad=[];
0069 showpcs(currpcs, mixedfilters, pixw, pixh)
0070 
0071 PCuse = setdiff(currpcs, PCbad);
0072 showpcs(PCuse, mixedfilters, pixw, pixh)
0073 
0074 fprintf('  Retaining PCs in the range [Klow - Khigh] = [%d - %d].\n', PCf,PCl)
0075 
0076 function showpcs(usepcs, Efull, pixw, pixh)
0077 
0078 if nargin<3
0079     fprintf('Assuming movie frames are square.\n')
0080     pixw = sqrt(size(Efull,1));
0081     pixh = sqrt(size(Efull,1));
0082 end
0083 if isempty(usepcs)
0084     usepcs = [1:size(Efull,2)];
0085 end
0086 
0087 if ndims(Efull)>=3
0088     Efull = reshape(Efull, pixw*pixh, []);
0089 end
0090 for j=usepcs
0091     Efull(:,j) = zscore(Efull(:,j));
0092 end
0093 pcs = reshape(Efull(:,usepcs), pixw, pixh, []);
0094 pcs = permute(pcs, [1, 2, 4, 3]);
0095 montage(pcs)
0096 colormap(hot)
0097 axis on
0098 xl = xlim;
0099 yl = ylim;
0100 nw = ceil(xl(2)/pixh)-1;
0101 nh = ceil(yl(2)/pixw)-1;
0102 set(gca,'YTick',[pixw:pixw:yl(2)],'YTickLabel',  num2str(usepcs(min([0:nh]*nw+1, length(usepcs)))'), ...
0103     'XTick',[pixh:pixh:xl(2)], ...
0104     'XTickLabel',num2str(usepcs([(nh-1)*nw+1:length(usepcs)])'), 'XAxisLocation','bottom','LineWidth',2)
0105 grid on
0106 formataxes
0107 caxis([-1,1]*7)
0108 
0109 function formataxes
0110 
0111 set(gca,'FontSize',12,'FontWeight','bold','FontName','Helvetica','LineWidth',2,'TickLength',[1,1]*.02,'tickdir','out')
0112 set(gcf,'Color','w','PaperPositionMode','auto')