Looking-to-cue-sides-and-modulation-of-value-in-OFC/code/support_functions/vecdiffsubsamp.m

function [vdg0nm, vdnmd, vdnmn, vd, vdg0, fpr, tpr, auc]=vecdiffsubsamp(x,y,nrecalls)
if nargin<3; nrecalls=1; end
if ~isvector(x) || ~isvector(y)
    error('inputs should be 1-D vectors');
end
x=rmnans(x(:));
y=rmnans(y(:));
Nx=length(x);
Ny=length(y);
vdg0nm=nan;
vdnmd=nan;
vdnmn=nan;
vd=nan;
vdg0=nan;
if Nx>0 && Ny>0
    zd=linspace(min([x(:); y(:)]), max([x(:); y(:)]), 1e3);
    if Nx==Ny
        xr=x;
        yr=y';
    elseif Ny>Nx
        nfolds=ceil(Ny/Nx);
        iy=mod(randperm(Nx*nfolds*nrecalls),Ny)+1;
        xr=repmat(x',Nx,1,nfolds*nrecalls);
        yr=reshape(repmat(y(iy)',1,Nx),Nx,Nx,nfolds*nrecalls);
    elseif Nx>Ny
        nfolds=ceil(Nx/Ny);
        ix=mod(randperm(Ny*nfolds*nrecalls),Nx)+1;
        xr=reshape(repmat(x(ix)',1,Ny),Ny,Ny,nfolds*nrecalls);
        yr=repmat(y',Ny,1,nfolds*nrecalls);
    end
    vd=xr-yr;
    [fpr,tpr,auc]=get_fpr_tpr_auc(x(:),y(:),zd);
    vd_nans=nan(size(vd));
    vd_nans(~isnan(vd))=1;
    vdg0=(vd>0).*vd_nans;
    vdg0nm=nanmean(vdg0,'all');
    vdnmd=nanmedian(vd,'all');
    vdnmn=nanmean(vd,'all');
end
end