isi

 Calculate the inter-spike-interval histogram                 
    Usage: [N,B,E] = isi(data,T,err,Nbins,plt)
    
 Input:                                                       
 Note that all times have to be consistent. 

 data   - structure array of spike times  (required)             
 T      - time interval of interest (default all)             
 err    - 0 for no error bars, 1 for jackknife errors
 
 Nbins  - number of bins in the isi                           
                                                              
 Output:                                                      
                                                              
 N      - count in bins                                       
 B      - bin centres                                            
 E      - errorbar (this is 2 sig deviation                   
          calculated using a jackknife over trials)

0001 function[N,B,E] = isi(data,T,err,Nbins,plt)
0002 % Calculate the inter-spike-interval histogram
0003 %    Usage: [N,B,E] = isi(data,T,err,Nbins,plt)
0004 %
0005 % Input:
0006 % Note that all times have to be consistent.
0007 %
0008 % data   - structure array of spike times  (required)
0009 % T      - time interval of interest (default all)
0010 % err    - 0 for no error bars, 1 for jackknife errors
0011 %
0012 % Nbins  - number of bins in the isi
0013 %
0014 % Output:
0015 %
0016 % N      - count in bins
0017 % B      - bin centres
0018 % E      - errorbar (this is 2 sig deviation
0019 %          calculated using a jackknife over trials)
0020 
0021 
0022 if nargin < 1; error('I need data!'); end
0023 data=padNaN(data); % create a zero padded data matrix from input structural array
0024 data=data'; % transposes data to get it in a form compatible with Murray's routine
0025 if nargin < 2; T = [min(data(:,1)) max(max(data))]; end
0026 if nargin < 3; err = 0;end
0027 if nargin < 4; Nbins = -1; end
0028 if nargin < 5; plt = 'r'; end
0029 
0030 if isempty(T); T = [min(min(data)) max(max(data))]; end
0031 if isempty(err); err = 0;end
0032 if isempty(Nbins); Nbins = -1; end
0033 if isempty(plt); plt = 'r'; end
0034 
0035 %  get the number of intervals in each trial and the indices of spike times
0036 %  that are kept
0037 
0038 NT = length(data(1,:)); % number of trials
0039 NI=zeros(1,NT);
0040 index(1:NT)=struct('keep',[]);
0041 for n=1:NT
0042   indx = find(data(:,n) >=  T(1) & data(:,n) <=  T(2) ... 
0043                                  & ~isnan(data(:,n)));
0044   if isempty(indx)
0045     NI(n) = 0;
0046   else
0047     NI(n) = length(indx)-1;
0048     index(n).keep=indx;
0049   end 
0050 end
0051 
0052 
0053 % calculate intervals...
0054 
0055 I = zeros(NT,max(NI));
0056 IT = [];
0057 for n=1:NT
0058   I(n,1:NI(n)) = diff(data(index(n).keep,n));
0059   IT = [IT I(n,1:NI(n))];
0060 end
0061 
0062 Mx = max(IT);
0063 if Nbins == -1
0064   Nbins = floor(sum(NI)/30);
0065   Med = median(IT);
0066   Nbins = max(floor(Nbins*Mx/Med),10);
0067 end
0068 
0069 B = linspace(0,Mx,Nbins);
0070 
0071 N = zeros(NT,Nbins);
0072 for n=1:NT
0073   N(n,:) = hist(I(n,1:NI(n)),B);
0074 end
0075 
0076 % answer...
0077 
0078 if NT > 1;Ns = sum(N)/NT;else Ns = N;end
0079 if ~strcmp(plt,'n')
0080   bar(B,NT*Ns);
0081 end
0082 
0083 % Jackknife iver trials to estimate std...
0084 
0085 if NT > 4 && err == 1
0086   MN = 0;
0087   SN = 0;
0088   for n=1:NT
0089     JK = (NT*Ns - N(n,:))/(NT-1);
0090     MN = MN + JK;
0091     SN = SN + JK.^2;   
0092   end  
0093   MN = MN/NT;
0094   SN = SN/NT;
0095   E = sqrt((NT-1)*(SN - MN.^2));
0096   if ~strcmp(plt,'n')
0097     hold on
0098     errorbar(B,NT*Ns,NT*2*E,'r-')
0099     hold off
0100   end
0101 end
0102 N = NT*Ns;
0103 
0104 
0105 
0106 
0107 
0108 
0109 
0110