Home > NoiseTools > nt_dprime_old.m

# nt_dprime_old

## PURPOSE

[d,err_rate,perm]=nt_dprime(x,y,jd_flag) - calculate d' (discriminability) of two distributions

## SYNOPSIS

function [d,err_rate,perm_rate]=nt_dprime(x,y,jd_flag)

## DESCRIPTION

    [d,err_rate,perm]=nt_dprime(x,y,jd_flag) - calculate d' (discriminability) of two distributions

      d: discriminablity index
      err_rate: error rate for linear discrimination
      perm_rate: rate for permuted data

      x, y: data (column vectors or matrices)
      jd_flag: apply JD first

     NoiseTools

## CROSS-REFERENCE INFORMATION

This function calls:

- nt_cov [c,tw]=nt_cov(x,shifts,w) - time shift covariance
- nt_demean [y,mn]=nt_demean(x,w) - remove weighted mean over cols
- nt_dprime [d,e]=nt_dprime(x,y,jd_flag) - calculate d' (discriminability) of two distributions
- nt_dss0 [todss,pwr1,pwr2]=nt_dss0(c0,c1,keep1,keep2) - dss from covariance
- nt_mmat y=nt_mmat(x,m) - matrix multiplication (with convolution)

This function is called by:

## SOURCE CODE

    0001 function [d,err_rate,perm_rate]=nt_dprime(x,y,jd_flag)
    0002 %[d,err_rate,perm]=nt_dprime(x,y,jd_flag) - calculate d' (discriminability) of two distributions
    0003 %
    0004 %  d: discriminablity index
    0005 %  err_rate: error rate for linear discrimination
    0006 %  perm_rate: rate for permuted data
    0007 %
    0008 %  x, y: data (column vectors or matrices)
    0009 %  jd_flag: apply JD first
    0010 %
    0011 % NoiseTools
    0012
    0013 NSTEPS=1000; % number of steps to find min of error rate
    0014 P=0.05; % theshold value for permutation test
    0015 NPERMUTE=1000; % number of trials for permutation test
    0016
    0017 if nargin<3; jd_flag=[]; end
    0018 if nargin<2; error('!'); end
    0019 if iscell(x)
    0020     xx=[]; yy=[];
    0021     for iCell=1:numel(x);
    0022         xx=[xx; x{iCell}];
    0023         yy=[yy; y{iCell}];
    0024     end
    0025     x=xx; y=yy;
    0026 end
    0027 if size(x,2) ~= size(y,2); error('!'); end
    0028
    0029 if jd_flag;
    0030     c0=nt_cov(nt_demean(x))+nt_cov(nt_demean(y));
    0031     c1=nt_cov(mean(x)-mean(y));
    0032     todss=nt_dss0(c0,c1);
    0033     x=nt_mmat(x,todss);
    0034     y=nt_mmat(y,todss);
    0035 end
    0036
    0037 d=abs(mean(x)-mean(y)) ./ sqrt((var(x)+var(y))/2);
    0038
    0039 % make sure that y>x
    0040 for iChan=1:size(x,2)
    0041     if mean(x(:,iChan))>mean(y(:,iChan));
    0042         x(:,iChan)=-x(:,iChan);
    0043         y(:,iChan)=-y(:,iChan);
    0044     end
    0045 end
    0046
    0047 if nargout>1; % error rate
    0048     err_rate=[];
    0049     for iChan=1:size(x,2)
    0050         min_error=1;
    0051         for thresh=linspace(mean(x(:,iChan)),mean(y(:,iChan)),NSTEPS);
    0052             x2y=sum(x(:,iChan)>thresh);
    0053             y2x=sum(y(:,iChan)<thresh);
    0054             nErr=x2y+y2x;
    0055             min_error=min(min_error,nErr/(size(x,1)+size(y,1)));
    0056         end
    0057         err_rate(iChan)=min_error;
    0058     end
    0059 end
    0060
    0061 if nargout>2; %permutation test
    0062     perm_rate=[];
    0063     for iChan=1:size(x,2)
    0064         for iPermute=1:NPERMUTE
    0065             if rem(iPermute,100)==1; disp(iPermute); end
    0066             % scramble between x and y:
    0067             z=[x(:,iChan);y(:,iChan)];
    0068             z=z(randperm(size(z,1)));
    0069             xx=z(1:size(x,1));
    0070             yy=z(size(x,1)+1:end);
    0071             min_error=1;
    0072             % scan criterion for minimum error
    0073             for thresh=linspace(mean(x(:,iChan)),mean(y(:,iChan)),NSTEPS);
    0074                 x2y=sum(xx>thresh);
    0075                 y2x=sum(yy<thresh);
    0076                 nErr=x2y+y2x;
    0077                 min_error=min(min_error,nErr/(size(x,1)+size(y,1)));
    0078             end
    0079             min_errors(iPermute)=min_error;
    0080         end
    0081         % find 5th percentile of distribution of error rates
    0082         min_errors=sort(min_errors);
    0083         min_errors=min_errors(1:round(NPERMUTE*P));
    0084         perm_rate(iChan)=min_errors(end);
    0085     end
    0086 end
    0087
    0088
    0089 % test code
    0090 if 0
    0091     x=randn(10000,1);
    0092     y=1+randn(10000,1);
    0093     figure(1); clf
    0094     t=-3:0.1:4;
    0095     plot(t,hist(x,t));
    0096     hold on;
    0097     plot(t,hist(y,t), 'r');
    0098     [d,e,p]=nt_dprime(x,y);
    0099     disp(['d'': ', num2str(d)]);
    0100     disp(['error: ', num2str(e)]);
    0101     disp(['95th percentile of error: ', num2str(p)]);
    0102 end