zoukankan      html  css  js  c++  java
  • perfcurve.m

    function [X,Y,T,auc,optrocpt,subY,subYnames] = ...
        perfcurve(labels,scores,posClass,varargin)
    %PERFCURVE Compute Receiver Operating Characteristic (ROC) curve or other
    %   performance curve for classifier output.
    %  
    %   [X,Y] = PERFCURVE(LABELS,SCORES,POSCLASS) computes a ROC curve for a
    %   vector of classifier predictions SCORES given true class labels,
    %   LABELS. The labels can be a numeric vector, logical vector, character
    %   matrix, cell array of strings or categorical vector (see help for
    %   groupingvariable). SCORES is a numeric vector of scores returned by a
    %   classifier for some data. This vector must have as many elements as
    %   LABELS does. POSCLASS is the positive class label (scalar). POSCLASS
    %   can be numeric (for numeric LABELS), logical (for logical LABELS), a
    %   character string (for character LABELS), a character string or a
    %   cellstr scalar (when LABELS are a cell array of strings), and a
    %   categorical scalar (for categorical LABELS). The specified positive
    %   class must be in the array of input labels. The returned values X and Y
    %   are coordinates for the performance curve and can be visualized  with
    %   PLOT(X,Y). By default, X is false positive rate, FPR, (equivalently,
    %   fallout, or 1-specificity) and Y is true positive rate, TPR,
    %   (equivalently, recall, or sensitivity).
    %  
    %   [X,Y,T] = PERFCURVE(LABELS,SCORES,POSCLASS) returns an array T of
    %   thresholds on classifier scores for the computed values of X and Y. It
    %   has the same number of rows as X and Y. For each threshold, TP is the
    %   count of true positive observations with scores greater or equal to
    %   this threshold, and FP is the count of false positive observations with
    %   scores greater or equal to this threshold. PERFCURVE defines negative
    %   counts, TN and FN, in a similar way and sorts the thresholds in the
    %   descending order which corresponds to the ascending order of positive
    %   counts. For the M distinct thresholds found in the array of scores,
    %   PERFCURVE returns the X, Y and T arrays with M+1 rows. PERFCURVE
    %   sets elements T(2:M+1) to the distinct thresholds, and T(1) replicates
    %   T(2). By convention, T(1) represents the highest 'reject all' threshold
    %   and PERFCURVE computes the corresponding values of X and Y for TP=0 and
    %   FP=0. T(end) is the lowest 'accept all' threshold for which TN=0 and
    %   FN=0.
    %  
    %   [X,Y,T,AUC] = PERFCURVE(LABELS,SCORES,POSCLASS) returns the area under
    %   curve (AUC) for the computed values of X and Y. Unless you specify
    %   XVALS or TVALS, PERFCURVE computes AUC using the returned X and Y
    %   values. If XVALS or TVALS is a numeric array, PERFCURVE computes AUC
    %   using X and Y values found from all distinct scores in the interval
    %   specified by the smallest and largest elements of XVALS or TVALS. For
    %   example, if XVALS is a numeric array, PERFCURVE finds X values for all
    %   distinct thresholds and uses a subset of these (with corresponding Y
    %   values) between MIN(XVALS) and MAX(XVALS) to compute AUC. The function
    %   uses trapezoidal approximation to estimate the area.
    %  
    %   If the first or last value of X or Y are NaN's, PERFCURVE removes them
    %   to allow calculation of AUC. This takes care of criteria that produce
    %   NaN's for the special 'reject all' or 'accept all' thresholds, for
    %   example, positive predictive value (PPV) or negative predictive value
    %   (NPV).
    %
    %   [X,Y,T,AUC] = PERFCURVE(LABELS,SCORES,POSCLASS) also returns pointwise
    %   confidence bounds for the computed values X, Y, T, and AUC if you
    %   supply cell arrays for LABELS and SCORES or set NBOOT to a positive
    %   integer. To compute the confidence bounds, PERFCURVE uses either
    %   vertical averaging (VA) or threshold averaging (TA). The returned
    %   values Y are an M-by-3 array in which the 1st element in every row
    %   gives the mean value, the 2nd element gives the lower bound and the 3rd
    %   element gives the upper bound. The returned AUC is a row-vector with 3
    %   elements following the same convention. For VA, the returned values T
    %   are an M-by-3 array and X is a column-vector. For TA, the returned
    %   values X are an M-by-3 matrix and T is a column-vector. 
    %
    %   PERFCURVE computes confidence bounds using either cross-validation or
    %   bootstrap. If you supply cell arrays for LABELS and SCORES, PERFCURVE
    %   uses cross-validation and treats elements in the cell arrays as
    %   cross-validation folds. LABELS can be a cell array of numeric vectors,
    %   logical vectors, character matrices, cell arrays of strings or
    %   categorical vectors. All elements in LABELS must have the same type.
    %   SCORES is a cell array of numeric vectors. The cell arrays for LABELS
    %   and SCORES must have the same number of elements, and the number of
    %   labels in cell K must be equal to the number of scores in cell K for
    %   any K in the range from 1 to the number of elements in SCORES.
    %
    %   If you set NBOOT to a positive integer, PERFCURVE generates NBOOT
    %   bootstrap replicas to compute pointwise confidence bounds. You cannot
    %   supply cell arrays for LABELS and SCORES and set NBOOT to a positive
    %   integer at the same time.
    %
    %   PERFCURVE returns pointwise confidence bounds. It does not return
    %   a simultaneous confidence band for the entire curve.
    %
    %   If you use 'XCrit' or 'YCrit' options described below to set the
    %   criterion for X or Y to an anonymous function, PERFCURVE can only
    %   compute confidence bounds by bootstrap.
    %    
    %   [X,Y,T,AUC,OPTROCPT] = PERFCURVE(LABELS,SCORES,POSCLASS) returns the
    %   optimal operating point of the ROC curve as an array of size 1-by-2
    %   with FPR and TPR values for the optimal ROC operating point. OPTROCPT
    %   is computed only for the standard ROC curve and set to NaN's otherwise.
    %   To obtain the optimal operating point for the ROC curve, PERFCURVE
    %   first finds the slope, S, using 
    %          S = (cost(P|N)-cost(N|N))/(cost(N|P)-cost(P|P)) * N/P
    %   where cost(I|J) is the cost of assigning an observation of class J to
    %   class I, and P=TP+FN and N=TN+FP are the total observation counts in
    %   the positive and negative class, respectively. PERFCURVE then finds the
    %   optimal operating point by moving the straight line with slope S from
    %   the upper left corner of the ROC plot (FPR=0,TPR=1) down and to the
    %   right until it intersects the ROC curve.
    %  
    %   [X,Y,T,AUC,OPTROCPT,SUBY] = PERFCURVE(LABELS,SCORES,POSCLASS) returns
    %   an array of Y values for negative subclasses. If you only specify one
    %   negative class, SUBY is identical to Y. Otherwise SUBY is a matrix of
    %   size M-by-K where M is the number of returned values for X and Y, and K
    %   is the number of negative classes. PERFCURVE computes Y values by
    %   summing counts over all negative classes. SUBY gives values of the Y
    %   criterion for each negative class separately. For each negative class,
    %   PERFCURVE places a new column in SUBY and fills it with Y values for TN
    %   and FP counted just for this class.
    %  
    %   [X,Y,T,AUC,OPTROCPT,SUBY,SUBYNAMES] = PERFCURVE(LABELS,SCORES,POSCLASS,'NegClass',NEGCLASS)
    %   returns a cell array of negative class names. If you provide an input
    %   array, NEGCLASS, of negative class names, PERFCURVE copies it into
    %   SUBYNAMES. If you do not provide NEGCLASS, PERFCURVE extracts SUBYNAMES
    %   from input labels. The order of SUBYNAMES is the same as the order of
    %   columns in SUBY, that is, SUBY(:,1) is for negative class SUBYNAMES{1}
    %   etc.
    %  
    %   [X,Y] = PERFCURVE(LABELS,SCORES,POSCLASS,'PARAM1',val1,'PARAM2',val2,...) 
    %   specifies optional parameter name/value pairs:
    %  
    %     'NegClass' - List of negative classes. Can be either a numeric array
    %                  or an array of chars or a cell array of strings. By
    %                  default, NegClass is set to 'all' and all classes found
    %                  in the input array of labels that are not the positive
    %                  class are considered negative. If NegClass is a subset
    %                  of the classes found in the input array of labels,
    %                  observations with labels that do not belong to either
    %                  positive or negative classes are discarded.
    %  
    %     'XCrit' - Criterion to compute for X. This criterion must be a
    %               monotone function of the positive class score. The
    %               following criteria are supported: 
    %         TP    - number of true positives
    %         FN    - number of false negatives
    %         FP    - number of false positives
    %         TN    - number of true negatives
    %         TP+FP - sum of TP and FP
    %         RPP   = (TP+FP)/(TP+FN+FP+TN) rate of positive predictions
    %         RNP   = (TN+FN)/(TP+FN+FP+TN) rate of negative predictions
    %         accu  = (TP+TN)/(TP+FN+FP+TN) accuracy
    %         TPR, sens, reca = TP/(TP+FN) true positive rate, sensitivity, recall
    %         FNR, miss       = FN/(TP+FN) false negative rate, miss
    %         FPR, fall       = FP/(TN+FP) false positive rate, fallout
    %         TNR, spec       = TN/(TN+FP) true negative rate, specificity
    %         PPV, prec = TP/(TP+FP) positive predictive value, precision
    %         NPV       = TN/(TN+FN) negative predictive value
    %         ecost=(TP*COST(P|P)+FN*COST(N|P)+FP*COST(P|N)+TN*COST(N|N))/(TP+FN+FP+TN)
    %              expected cost
    %         In addition, you can define an arbitrary criterion by supplying
    %         an anonymous function of 3 arguments, (C,scale,cost), where C is
    %         a 2-by-2 confusion matrix, scale is a 2-by-1 array of class
    %         scales, and cost is a 2-by-2 misclassification cost matrix. See
    %         doc for Performance Curves for more info. Warning: some of these
    %         criteria return NaN values at one of the two special thresholds,
    %         'reject all' and 'accept all'.
    %   
    %     'YCrit' - Criterion to compute for Y. The same criteria as for X
    %               are supported. This criterion does not have to be a
    %               monotone function of the positive class score.
    %   
    %     'XVals' - Values for the X criterion. By default, XVals is unset
    %               and PERFCURVE computes X, Y and T values for all scores.
    %               You can set XVals to either 'all' or a numeric array. If
    %               XVals is set to 'all' and TVals is unset, PERFCURVE returns
    %               X, Y and T values for all scores and computes pointwise
    %               confidence bounds for Y and T using vertical averaging. If
    %               XVals is set to a numeric array, PERFCURVE returns X, Y and
    %               T values for the specified XVals and computes pointwise
    %               confidence bounds for Y and T at these XVals using vertical
    %               averaging. You cannot set XVals and TVals at the same time.
    %
    %     'TVals' - Thresholds for the positive class score. By default, TVals
    %               is unset and PERFCURVE computes X, Y and T values for all
    %               scores. You can set TVals to either 'all' or a numeric
    %               array. If TVals is set to 'all' or unset and XVals is
    %               unset, PERFCURVE returns X, Y and T values for all scores
    %               and computes pointwise confidence bounds for Y and X using
    %               threshold averaging. If TVals is set to a numeric array,
    %               PERFCURVE returns X, Y and T values for the specified
    %               thresholds and computes pointwise confidence bounds for Y
    %               and X at these thresholds using threshold averaging. You
    %               cannot set XVals and TVals at the same time. 
    %   
    %     'UseNearest' - 'on' to use nearest values found in the data instead
    %                    of the specified numeric XVals or TVals and 'off'
    %                    otherwise. If you specify numeric XVals and set
    %                    UseNearest to 'on', PERFCURVE returns nearest unique
    %                    values X found in the data, as well as corresponding
    %                    values of Y and T. If you specify numeric XVals and
    %                    set UseNearest to 'off', PERFCURVE returns these XVals
    %                    sorted. By default this parameter is set to 'on'. If
    %                    you compute confidence bounds by cross-validation or
    %                    bootstrap, this parameter is always 'off'.
    %
    %     'ProcessNaN' - This argument specifies how PERFCURVE processes NaN
    %                    scores. By default, it is set to 'ignore' and
    %                    observations with NaN scores are removed from the
    %                    data. If the parameter is set to 'addtofalse',
    %                    PERFCURVE adds observations with NaN scores to false
    %                    classification counts in the respective class. That
    %                    is, observations from the positive class are always
    %                    counted as false negative (FN), and observations from
    %                    the negative class are always counted as false
    %                    positive (FP). 
    %   
    %     'Prior' - Either string or array with 2 elements. It represents prior
    %               probabilities for the positive and negative class,
    %               respectively. Default is 'empirical', that is, prior
    %               probabilities are derived from class frequencies. If set to
    %               'uniform', all prior probabilities are set equal.
    %   
    %     'Cost'  - A 2-by-2 matrix of misclassification costs 
    %                   [C(P|P) C(N|P); C(P|N) C(N|N)] 
    %               where C(I|J) is the cost of misclassifying
    %               class J as class I. By default set to [0 0.5; 0.5 0].
    %
    %     'Alpha' - A numeric value between 0 and 1. PERFCURVE returns
    %               100*(1-Alpha) percent pointwise confidence bounds for X, Y,
    %               T and AUC. By default set to 0.05 for 95% confidence
    %               interval.
    %
    %     'Weights' - A numeric vector of non-negative observation weights.
    %                 This vector must have as many elements as SCORES or
    %                 LABELS do. If you supply cell arrays for SCORES and
    %                 LABELS and you need to supply WEIGHTS, you must supply
    %                 them as a cell array too. In this case, every element in
    %                 WEIGHTS must be a numeric vector with as many elements as
    %                 the corresponding element in SCORES:
    %                 NUMEL(WEIGHTS{1})==NUMEL(SCORES{1}) etc. To compute X, Y
    %                 and T or to compute confidence bounds by
    %                 cross-validation, PERFCURVE uses these observation
    %                 weights instead of observation counts. To compute
    %                 confidence bounds by bootstrap, PERFCURVE samples N out
    %                 of N with replacement using these weights as multinomial
    %                 sampling probabilities.
    %
    %     'NBoot' - Number of bootstrap replicas for computation of confidence
    %               bounds. Must be a positive integer. By default this
    %               parameter is set to zero, and bootstrap confidence bounds
    %               are not computed. If you supply cell arrays for LABELS and
    %               SCORES, this parameter must be set to zero because
    %               PERFCURVE cannot use both cross-validation and bootstrap to
    %               compute confidence bounds.
    %
    %     'BootType' - Confidence interval type used by BOOTCI to compute
    %                  confidence  bounds. You can specify any type supported
    %                  by BOOTCI. 'doc bootci' for more info. By default set to
    %                  'bca'.
    %
    %     'BootArg' - Optional input arguments for BOOTCI used to compute
    %                 confidence bounds. You can specify all arguments
    %                 supported by BOOTCI. 'doc bootci' for more info. Empty by
    %                 default.
    %     
    %   Example: Plot ROC curve for classification by logistic regression
    %      load fisheriris
    %      x = meas(51:end,1:2);        % iris data, 2 classes and 2 features
    %      y = (1:100)'>50;             % versicolor=0, virginica=1
    %      b = glmfit(x,y,'binomial');  % logistic regression
    %      p = glmval(b,x,'logit');     % get fitted probabilities for scores
    % 
    %      [X,Y] = perfcurve(species(51:end,:),p,'virginica');
    %      plot(X,Y)
    %      xlabel('False positive rate'); ylabel('True positive rate')
    %      title('ROC for classification by logistic regression')
    %
    %      % Obtain errors on TPR by vertical averaging
    %      [X,Y] = perfcurve(species(51:end,:),p,'virginica','nboot',1000,'xvals','all');
    %      errorbar(X,Y(:,1),Y(:,2)-Y(:,1),Y(:,3)-Y(:,1)); % plot errors
    %
    %   See also fitglm, fitcdiscr, fitNaiveBayes, fitctree, fitcknn, fitcsvm,
    %   fitensemble, TreeBagger, groupingvariable, bootci.
    
    %   Copyright 2008-2013 The MathWorks, Inc.
    
    
    % Prepare input parser
    p = inputParser;
    p.addRequired('labels', ...
        @(x) ~isempty(x) && (ischar(x) || isnumeric(x) || islogical(x) ...
        || iscell(x) || isa(x,'categorical')));
    p.addRequired('scores',@(x) ~isempty(x) && (isnumeric(x) || iscell(x)));
    p.addRequired('PosClass', ...
        @(x) ~isempty(x) && (ischar(x) || isnumeric(x) || islogical(x)...
        || iscellstr(x) || isa(x,'categorical')));
    p.addParamValue('NegClass','all', ...
        @(x) ~isempty(x) && (ischar(x) || isnumeric(x) || iscell(x)));    
    p.addParamValue('XCrit','FPR',@(x) ischar(x) || isa(x,'function_handle'));
    p.addParamValue('YCrit','TPR',@(x) ischar(x) || isa(x,'function_handle'));
    p.addParamValue('XVals','', ... 
        @(x) (ischar(x) && (strcmpi(x,'all') || isempty(x))) || ...
        (~isempty(x) && isnumeric(x)) );
    p.addParamValue('TVals','', ... 
        @(x) (ischar(x) && (strcmpi(x,'all') || isempty(x))) || ...
        (~isempty(x) && isnumeric(x)) );
    p.addParamValue('UseNearest','on',@(x) strcmpi(x,'on') || strcmpi(x,'off'));
    p.addParamValue('ProcessNaN','ignore', ...
        @(x) strcmpi(x,'ignore') || strcmpi(x,'addtofalse'));
    p.addParamValue('Prior','empirical', ...
        @(x) (ischar(x) && (strcmpi(x,'empirical') || strcmpi(x,'uniform'))) ...
        || (isnumeric(x) && numel(x)==2) );
    p.addParamValue('Cost',[0 0.5; 0.5 0], ...
        @(x) isnumeric(x) && size(x,1)==2 && size(x,2)==2);
    p.addParamValue('Weights',[],@(x) isempty(x) || isfloat(x) || iscell(x));
    p.addParamValue('NBoot',0,@(x) isnumeric(x) && isscalar(x) && x>=0);
    p.addParamValue('BootType','bca',@(x) ischar(x));
    p.addParamValue('BootArg',{},@(x) iscell(x));
    p.addParamValue('Alpha',0.05,@(x) isnumeric(x) && isscalar(x) && x>0 && x<1);
    p.FunctionName = 'perfcurve';
    
    % Parse inputs
    p.parse(labels,scores,posClass,varargin{:});
    negClass       = p.Results.NegClass;
    xCrit          = p.Results.XCrit;
    yCrit          = p.Results.YCrit;
    xVals          = p.Results.XVals;
    tVals          = p.Results.TVals;
    toNearest      = strcmpi(p.Results.UseNearest,'on');
    processNaN     = p.Results.ProcessNaN;
    prior          = p.Results.Prior;
    cost           = p.Results.Cost;
    weights        = p.Results.Weights;
    nboot          = ceil(p.Results.NBoot);
    boottype       = p.Results.BootType;
    bootarg        = p.Results.BootArg;
    alpha          = p.Results.Alpha;
    
    % Prepare data
    [scores,labels,weights,ncv] = preparedata(scores,labels,weights);
    
    % Check compatibility of thresholds and xvals arguments
    % By default use supplied thresholds for computing the curve
    useTVals = true; 
    if ~isempty(tVals) && ~isempty(xVals)
        error(message('stats:perfcurve:BothTandXsupplied'));
    end
    % If X values are supplied, use them to compute the curve
    if ~isempty(xVals)
        useTVals = false;
    end
    
    % Check if both cross-validation and bootstrap are requested
    docv = false;
    doboot = false;
    if ncv>0
        docv = true;
    end
    if nboot>0
        doboot = true;
    end
    if docv && doboot
        error(message('stats:perfcurve:BothCVandBootstrapRequested'));
    end
    nsub = max(ncv,nboot);
    
    % Stack all CV values into one long vector
    % ncvel is the number of elements per CV fold
    if docv
        [ncvel,scores,labels,weights] = stackcv(scores,labels,weights);
        ncvcum = [0; cumsum(ncvel)];
    end
    
    % Convert class labels to a cat array
    labels = nominal(labels);
    trueNames = categories(labels)'; % as a row
    if length(trueNames) < 2
        error(message('stats:perfcurve:NotEnoughClasses'));
    end
    
    % Check costs
    if (cost(2,1)-cost(2,2))<=0 || (cost(1,2)-cost(1,1))<=0
        error(message('stats:perfcurve:InvalidCost'));
    end
    
    % Sort scores in the descending order
    [sScores,sorted] = sort(scores,1,'descend');
    
    % Get class membership for instances:
    % W(i,j) is weight of observation i if observation i is from class j and 0
    % otherwise. 
    % Also, get negative class names.
    % W has the size of NxK,
    %   where N is the number of instances and K is the number of classes.
    % subYnames is a cell array of length K-1 with names of negative classes.
    % Column W(:,j) is for class subYnames{j-1}   (j>1)
    [W,subYnames] = membership(labels(sorted),weights(sorted),...
        posClass,negClass,trueNames);
    
    % Make Wcum, a matrix of cumulative weights in each class.
    % Adjust Wcum and scores using the specified behavior for NaN scores.
    % Output sorted distinct scores into sScores and corresponding rows into Wcum.
    % Wcum and output sScores do not have the same size as W and input sScores.
    % To access the full vector of scores, use scores(sorted).
    [Wcum,sScores] = makeccum(W,sScores,processNaN);
    
    % Get cumulative counts and sorted scores for each subset
    if docv
        subscores = cell(nsub,1);
        Wcumsub = cell(nsub,1);
        for isub=1:nsub
            tf = sorted>ncvcum(isub) & sorted<=ncvcum(isub+1);
            [Wcumsub{isub},subscores{isub}] = ...
                makeccum(W(tf,:),scores(sorted(tf)),processNaN);
        end
    end
    
    % Check that confidence bound computation by CV is not requested for
    % user-defined criteria
    if (isa(xCrit,'function_handle') || isa(yCrit,'function_handle')) && docv
        error(message('stats:perfcurve:UserCritConfBounds'));
    end
    
    % Determine criteria to compute
    % The 1st output is the function for computing criterion itself.
    % The 2nd output is the function for computing weight for this criterion.
    [fx,fwx] = makeCrit(xCrit);
    [fy,fwy] = makeCrit(yCrit);
    
    % Convert class probabilities into class scales.
    scale = classscale(Wcum,prior);
    
    % Define arrayfuns
    afx = @(tp,fn,fp,tn) fx([tp fn; fp tn],scale,cost);
    afy = @(tp,fn,fp,tn) fy([tp fn; fp tn],scale,cost);
    afwx = @(tp,fn,fp,tn) fwx([tp fn; fp tn],scale);
    afwy = @(tp,fn,fp,tn) fwy([tp fn; fp tn],scale);
    
    % Compute threshold indices
    uniqdiv = ~docv && ~doboot && toNearest;
    Ndiv = size(Wcum,1) - 1;
    if ischar(xVals) && ischar(tVals)
        div = (1:Ndiv)'; % Use all thresholds
    else
        if useTVals % Use user-defined thresholds
            [div,tVals] = tdiv(tVals,sScores,uniqdiv);
        else % Use user-defined X values
            [div,xVals] = xdiv(xVals,afx,Wcum,uniqdiv);
        end
    end
    
    % Compute the actual values for the specified criterion,
    %   (to be plotted on X axis),
    %   and associated TP and FP counts.
    [X,~,tpX,fpX] = Xvalues(div,Wcum,afx,[]);
    % If xVals are supplied by the user and if they do not need to be set to
    % nearest found values, use these xVals
    if ~uniqdiv && isnumeric(xVals)
        X = xVals(:);
    end
    
    % Compute criterion values associated with these thresholds 
    %   (to be plotted on Y axis)
    Y = Yvalues(tpX,fpX,Wcum,afy,[]);
    
    % Check if any criteria are NaN's besides those at special 'reject all' and
    % 'accept all' thresholds 
    special = (div==1 | div==Ndiv);
    if any(isnan(X(~special)))
        error(message('stats:perfcurve:BadXCritValue'));
    end
    if any(isnan(Y(~special)))
        error(message('stats:perfcurve:BadYCritValue'));
    end
    
    % Get thresholds from indices
    if nargout>2
        % If tVals are supplied by the user and if they do not need to be set
        % to nearest found values, use these tVals
        if ~uniqdiv && isnumeric(tVals)
            T = tVals(:);
        else
            T = thresholds(div,sScores,false);
        end
    end
    
    % Find numeric divisions for VA or TA
    tValsFixed = [];
    xValsFixed = [];
    if docv || doboot
        if useTVals
            if ischar(tVals)
                tValsFixed = thresholds(div,sScores,true);
            else
                tValsFixed = tVals;
            end
        else
            if ischar(xVals)
                xValsFixed = X;
            else
                xValsFixed = xVals;
            end
        end
    end
    
    % Compute confidence intervals
    ciX = [];
    ciY = [];
    ciT = [];
    if docv || doboot
        notifyClassAbsent(true); % Reset notification flag.
        if     docv % cross-validation
            % For CV, compute X, Y and T for supplied folds with weights
            [Xsub,wXsub,Ysub,wYsub,Tsub,wTsub] = ...
                xyt(xValsFixed,tValsFixed,subscores,Wcumsub,afx,afwx,afy,afwy);
            
            % X errors for TA
            if useTVals
                ciX = cvci(X,Xsub,wXsub,alpha);
            end
            
            % Score threshold errors for VA
            computeT = nargout>2 && ~useTVals;
            if computeT
                ciT = cvci(T,Tsub,wTsub,alpha);
            end
            
            % Always need Y errors
            ciY = cvci(Y,Ysub,wYsub,alpha);
        elseif doboot
            % Use bootci to compute confidence intervals. Use weights as
            % multinomial probabilities for bootstrap replica generation and
            % use observation counts (not observation weights!) to compute
            % performance curves on these replicas.
            bootfun = @(idx) ...
                oneBootXYT(idx,scores(sorted),W>0,xValsFixed,tValsFixed,processNaN,afx,afy);
            ci = bootci(nboot,{bootfun (1:size(W,1))'},...
                'weights',sum(W,2),'alpha',alpha,'type',boottype,bootarg{:});        
            
            % X and Y errors for TA
            if useTVals
                if length(tValsFixed)>1
                    ciX = [ci(1,:,1)' ci(2,:,1)'];
                    ciY = [ci(1,:,2)' ci(2,:,2)'];
                else
                    ciX = ci(:,1)';
                    ciY = ci(:,2)';
                end
            % Y errors for VA
            else
                if length(xValsFixed)>1
                    ciY = [ci(1,:,2)' ci(2,:,2)'];
                else
                    ciY = ci(:,2)';
                end
            end
            
            % Score threshold errors for VA
            computeT = nargout>2 && ~useTVals;
            if computeT
                if length(xValsFixed)>1
                    ciT = [ci(1,:,1)' ci(2,:,1)'];
                else
                    ciT = ci(:,1)';
                end
            end
        end
    end
    
    % Compute area under curve
    if nargout>3
        % Set range and either 'x' or 't' to specify what range
        XorTrange = [];
        mode = '';
        if     ~ischar(xVals)
            mode = 'x';
            XorTrange = sort([xVals(1) xVals(end)]);
        elseif ~ischar(tVals)
            mode = 't';
            XorTrange = sort([tVals(1) tVals(end)]);
        end
    
        % Get area
        auc = AUC(XorTrange,mode,sScores,Wcum,afx,afy);
    
        % Get error
        if docv || doboot
            if     docv
                subauc = zeros(1,nsub);
                Wtotsub = zeros(1,nsub);
                for isub=1:nsub
                    [subauc(isub),Wtotsub(isub)] = ...
                        AUC(XorTrange,mode,subscores{isub},Wcumsub{isub},afx,afy);
                end
                ciauc = cvci(auc,subauc,Wtotsub,alpha);
            elseif doboot
                bootfun = @(idx) ...
                    oneBootAUC(idx,scores(sorted),W>0,XorTrange,mode,processNaN,afx,afy);
                ciauc = bootci(nboot,{bootfun (1:size(W,1))'},...
                    'weights',sum(W,2),'alpha',alpha,'type',boottype,bootarg{:});
                ciauc = ciauc';
            end
            auc = [auc ciauc];
        end
    end
    
    % Find optimal operation point for the standard ROC curve
    if nargout>4
        isroc = (strcmpi(xCrit,'FPR') || strcmpi(xCrit,'fall')) && ...
            (strcmpi(yCrit,'TPR') || strcmpi(yCrit,'sens') || strcmpi(yCrit,'reca'));
        if isroc
            optrocpt = findoptroc(X,Y,Wcum,scale,cost);
        else % Not a standard ROC curve.
            optrocpt = NaN(1,2);
        end
    end
    
    % Compute criterion values for individual negative classes
    if nargout>5
        subY = subYvalues(tpX,fpX,div,Wcum,afy);
    end
    
    % Include confidence intervals if they were computed
    if ~isempty(ciX)
        X = [X ciX];
    end
    if ~isempty(ciY)
        Y = [Y ciY];
    end
    if ~isempty(ciT)
        T = [T ciT];
    end
    end
    
    
    function [W,negClassNames] = membership(sLabels,sWeights,posClass,negClass,trueNames)
    % Find the positive class. Must have exactly one.
    posClass = cellstr(nominal(posClass));
    if length(posClass)>1
        error(message('stats:perfcurve:TooManyPositiveClasses'));
    end
    if ~ismember(posClass,trueNames)
        error(message('stats:perfcurve:PositiveClassNotFound'));
    end
    
    % Check negative class labels
    if strcmpi(negClass,'all')
        negClass = nominal(trueNames);
        [~,posClassLoc] = ismember(posClass,cellstr(negClass));
        negClass(posClassLoc) = [];
        negClass = nominal(negClass);
    else
        negClass = nominal(negClass);
        tf = ismember(cellstr(negClass),trueNames);
        if any(~tf)
            error(message('stats:perfcurve:NegativeClassNotFound'));
        end
        tf = ismember(posClass,cellstr(negClass));
        if tf
            error(message('stats:perfcurve:PositiveAndNegativeClassesOverlap'));
        end
    end
    nNeg = length(negClass);
    
    % Names of selected negative classes
    negClassNames = categories(negClass)'; % as a row
    
    % Check for duplicate entries
    if nNeg~=length(negClassNames)
        error(message('stats:perfcurve:DuplicateNegativeClasses'));
    end
    
    % Fill out the membership matrix
    % 1st column is for the positive class.
    % Columns 2:end are for negative classes.
    negClass = cellstr(negClass);
    C = false(length(sLabels),1+nNeg);
    C(:,1) = ismember(sLabels,posClass);
    for i=1:nNeg
        C(:,i+1) = ismember(sLabels,negClass(i));
    end
    
    % Get weighted membership matrix
    W = bsxfun(@times,C,sWeights);
    end
    
    
    function [Wcum,scores] = makeccum(W,scores,processNaN)
    % Discard instances that do not belong to any class
    idxNone = ~any(W,2);
    W(idxNone,:) = [];
    scores(idxNone) = [];
    
    % Get rid of NaN's in scores
    Wnanrow = zeros(1,size(W,2));
    idxNaN = isnan(scores);
    if strcmpi(processNaN,'addtofalse')
        if ~isempty(idxNaN)
            Wnanrow = sum(W(idxNaN,:),1);
        end
    end
    scores(idxNaN) = [];
    W(idxNaN,:) = [];
    
    % Make a matrix of counts with NaN instances included
    Wnan = zeros(size(W,1)+2,size(W,2));
    Wnan(1,2:end) = Wnanrow(2:end);% FP (always accepted)
    Wnan(2:end-1,:) = W;
    Wnan(end,1) = Wnanrow(1);% FN (always rejected)
    
    % Compute cumulative counts in each class
    Wcum = cumsum(Wnan,1);
    
    % Compact Wcum in case of identical scores
    idxEq = find( scores(1:end-1) < scores(2:end) + ...
        max([eps(scores(1:end-1)) eps(scores(2:end))],[],2) );
    Wcum(idxEq+1,:) = [];
    scores(idxEq) = [];
    end
    
    
    function scale = classscale(Wcum,prior)
    scale = zeros(2,1);
    Wpos = Wcum(end,1);
    Wneg = sum(Wcum(end,2:end),2);
    if ischar(prior) && strcmpi(prior,'empirical')
        scale = ones(2,1);
    end
    if ischar(prior) && strcmpi(prior,'uniform')
        prior = ones(2,1);
    end
    if isnumeric(prior)
        if any(prior<=0)
            error(message('stats:perfcurve:NonPositivePriors'));
        end
        scale(1) = prior(1)*Wneg;
        scale(2) = prior(2)*Wpos;
        scale = scale/sum(scale);
    end
    end
    
    
    function [f,wf] = makeCrit(crit)
    % If this is a user-supplied function, just return it
    if isa(crit,'function_handle')
        f = crit;
        wf = @(C,scale) sum(scale.*sum(C,2));
        return;
    end
    
    % Make the function given criterion name
    switch lower(crit)
        case 'tp'
            f = @(C,scale,cost) scale(1)*C(1,1);
            wf = @(C,scale) 1;
        case 'fn'
            f = @(C,scale,cost) scale(1)*C(1,2);
            wf = @(C,scale) 1;
        case 'fp'
            f = @(C,scale,cost) scale(2)*C(2,1);
            wf = @(C,scale) 1;
        case 'tn'
            f = @(C,scale,cost) scale(2)*C(2,2);
            wf = @(C,scale) 1;
        case 'tp+fp'
            f = @(C,scale,cost) sum(scale.*C(:,1),1);
            wf = @(C,scale) 1;
        case 'rpp'
            f = @(C,scale,cost) sum(scale.*C(:,1),1) / sum(scale.*sum(C,2));
            wf = @(C,scale) sum(scale.*sum(C,2));
        case 'rnp'
            f = @(C,scale,cost) sum(scale.*C(:,2),1) / sum(scale.*sum(C,2));
            wf = @(C,scale) sum(scale.*sum(C,2));
        case 'accu'
            f = @(C,scale,cost) ...
                (scale(1)*C(1,1)+scale(2)*C(2,2)) / sum(scale.*sum(C,2));
            wf = @(C,scale) sum(scale.*sum(C,2));
        case {'tpr','sens','reca'}
            f = @(C,scale,cost) C(1,1) / sum(C(1,:),2);
            wf = @(C,scale) sum(C(1,:),2);
        case {'fnr','miss'}
            f = @(C,scale,cost) C(1,2) / sum(C(1,:),2);
            wf = @(C,scale) sum(C(1,:),2);
        case {'fpr','fall'}
            f = @(C,scale,cost) C(2,1) / sum(C(2,:),2);
            wf = @(C,scale) sum(C(2,:),2);
        case {'tnr','spec'}
            f = @(C,scale,cost) C(2,2) / sum(C(2,:),2);
            wf = @(C,scale) sum(C(2,:),2);
        case {'ppv','prec'}
            f = @(C,scale,cost) scale(1)*C(1,1) / sum(scale.*C(:,1));
            wf = @(C,scale) sum(scale.*C(:,1));
        case 'npv'
            f = @(C,scale,cost) scale(2)*C(2,2) / sum(scale.*C(:,2));
            wf = @(C,scale) sum(scale.*C(:,2));
        case 'ecost'
            f = @(C,scale,cost) ...
                sum(scale.*sum(C.*cost,2)) / sum(scale.*sum(C,2));
            wf = @(C,scale) sum(scale.*sum(C,2));
        otherwise
            error(message('stats:perfcurve:UnknownXYCriterion'));
    end
    end
    
    
    function inrange = applyrange(X,T,XorTrange,mode)
    inrange = true(length(X),1);
    if ~isempty(XorTrange)
        if numel(XorTrange)~=2
            error(message('stats:perfcurve:InvalidXorTRange'));
        end
        XorTrange = sort(XorTrange);
        if     mode=='x'
            inrange = (X>=XorTrange(1) & X<=XorTrange(2));
        elseif mode=='t'
            inrange = (T>=XorTrange(1) & T<=XorTrange(2));
        end
        if isempty(inrange)
            error(message('stats:perfcurve:XorTRangeTooRestrictive'));
        end
    end
    end
    
    
    function T = thresholds(div,scores,shiftRejectAll)
    % Init
    T = zeros(length(div),1);
    
    % First threshold
    isone = div==1;
    T(isone) = scores(1);
    if shiftRejectAll
        T(isone) = T(isone) + eps(scores(1));
    end
    
    % Normal thresholds
    T(~isone) = scores(div(~isone)-1);
    end
    
    
    function increasing = monotone(vals)
    % Allow only monotone criteria.
    % By default, assume a criterion that monotonously increases as
    %   the predicted score in the positive class decreases. Otherwise, swap.
    % 'increasing' is a flag that shows in what order values are sorted.
    vals = vals(~isnan(vals));
    increasing = 1;
    if isempty(vals)
        return;
    end
    if any(vals(1:end-1)>vals(2:end))
        increasing = -1;
        if any(vals(1:end-1)<vals(2:end))
            error(message('stats:perfcurve:NonMonotoneXCriterion'));
        end
    end
    end
    
    %{
    % This routine does the same thing as the uncommented one below. It is
    % simpler but could be much slower for huge allVals.
    function div = finddiv(divVals,allVals,increasing)
    Nall = length(allVals);
    Nval = length(divVals);
    div = ones(Nval,1);
    for i=1:Nall
        div(increasing*divVals >= increasing*allVals(i)) = i;
    end
    end
    %}
    
    % Find allVals indices for values supplied in divVals.
    % allVals is a sorted array of all know values, and divVals is a sorted
    % array of values for which we want to find indices in allVals.
    % increasing is the sort order: +1 for ascending and -1 for descending.
    % Assume that divVals is not too large and allVals can be huge.
    % This routine works if allVals is a subset of divVals as well (case
    % relevant for cross-validation or bootstrap).
    function div = finddiv(divVals,allVals,increasing)
    % Init
    Ndiv = length(divVals);
    div = ones(Ndiv,1);
    
    % Set current search starting index
    iAll = 1;
    
    % Main loop
    for iDiv=1:Ndiv
        % Find index in the list of all known values
        thisAll = find(increasing*allVals(iAll:end) <= increasing*divVals(iDiv),1,'last');
        if isempty(thisAll)
            continue;
        end
        
        % Update starting index for the list of all known values
        iAll = iAll + thisAll - 1;
    
        % Fill found divisions
        div(iDiv:end) = iAll;
    end
    end
    
    
    % Find division indices along X for supplied xVals
    function [divX,xVals] = xdiv(xVals,afx,Wcum,uniqdiv)
    % Get counts for positive and negative classes
    Nrow = size(Wcum,1) - 1;
    Pcum = Wcum(:,1);
    Ncum = sum(Wcum(:,2:end),2);
    wP = Pcum(end);
    wN = Ncum(end);
    
    % Get all possible values of the criterion
    allVals = arrayfun(afx,Pcum(1:Nrow),wP-Pcum(1:Nrow),Ncum(1:Nrow),wN-Ncum(1:Nrow));
    
    % Do criterion values increase or decrease vs predicted scores?
    increasing = monotone(allVals);
    
    % Sort input values
    xVals = increasing*sort(increasing*xVals);
    
    % Find indices of thresholds.
    divX = finddiv(xVals,allVals,increasing);
    
    % Make indices unique if necessary 
    % and include the first threshold below 'accept all'.
    if uniqdiv
        if increasing*allVals(1)>=increasing*xVals(1)
            divX = [1; divX];
        end
        divX = unique(divX);
    end
    end
    
    
    % Find division indices along positive class score for supplied thresholds
    function [divT,tVals] = tdiv(tVals,scores,uniqdiv)
    % Sort thresholds
    tVals = sort(tVals,'descend');
    
    % Find indices of thresholds.
    % Scores have been sorted in the descending order.
    divT = finddiv(tVals,scores,-1);
    
    % Offset indices of non-special thresholds by 1 to account for the NaN row
    % on top of Ccum. 
    divT = divT + 1;
    divT(tVals>scores(1)) = 1;
    
    % Make indices unique if necessary
    if uniqdiv
        divT = unique(divT);
    end
    end
    
    
    % Get X values for given division indices
    function [valX,wX,tpX,fpX] = Xvalues(div,Wcum,afx,afwx)
    % Get counts for positive and negative classes
    Pcum = Wcum(:,1);
    Ncum = sum(Wcum(:,2:end),2);
    wP = Pcum(end);
    wN = Ncum(end);
    Nrow = size(Pcum,1)-1;
    
    % Check that division index does not exceed array size
    if any(div>Nrow)
        error(message('stats:perfcurve:BadDivisionIndex'));
    end
    
    % Get TP and FP counts for chosen threshold indices
    tpX = Pcum(div);
    fpX = Ncum(div);
    
    % valX is the corresponding array of criterion values
    valX = arrayfun(afx,tpX,wP-tpX,fpX,wN-fpX);
    
    % wX is the corresponding array of weights
    wX = [];
    if ~isempty(afwx)
        wX = arrayfun(afwx,tpX,wP-tpX,fpX,wN-fpX);
    end
    end
    
    
    function [valY,wY] = Yvalues(tpX,fpX,Wcum,afy,afwy)
    % Get number of instances in the positive and negative classes
    wP = Wcum(end,1);
    wN = sum(Wcum(end,2:end),2);
    
    % Compute Y criterion for the total of all negative classes
    valY = arrayfun(afy,tpX,wP-tpX,fpX,wN-fpX);
    
    % Compute corresponding weights
    wY = [];
    if ~isempty(afwy)
        wY = arrayfun(afwy,tpX,wP-tpX,fpX,wN-fpX);
    end
    end
    
    
    function subY = subYvalues(tpX,fpX,divX,Wcum,afy)
    % If only one negative class, it is accounted for by Yvalues
    nNegClass = size(Wcum,2)-1;
    if nNegClass < 2
        subY = Yvalues(tpX,fpX,Wcum,afy,[]);
        return;
    end
    
    % Compute Y criteria for negative classes separately
    wP = Wcum(end,1);
    subY = zeros(length(tpX),nNegClass);
    for i=1:nNegClass
        wN = Wcum(end,i+1);
        if wN==0
            error(message('stats:perfcurve:BadClassCounts', i));
        end
        fpX = Wcum(divX,i+1);
        subY(:,i) = arrayfun(afy,tpX,wP-tpX,fpX,wN-fpX);
    end
    end
    
    
    % Compute area under curve and associated weight
    function [auc,wauc] = AUC(XorTrange,mode,sScores,Wcum,afx,afy)
    % If not all thresholds have been found, find all and apply range.
    Ndiv = size(Wcum,1)-1;
    div = (1:Ndiv)';
    [X,~,tpx,fpx] = Xvalues(div,Wcum,afx,[]);
    Y = Yvalues(tpx,fpx,Wcum,afy,[]);
    T = thresholds(div,sScores,false);
    
    % Apply range
    inrange = applyrange(X,T,XorTrange,mode);
    needTrimFirst = true;
    if ~inrange(1)
        needTrimFirst = false;
    end
    needTrimLast = true;
    if ~inrange(end)
        needTrimLast = false;
    end
    X = X(inrange);
    Y = Y(inrange);
    icum = find(inrange,1,'last');
    if isempty(icum)
        wauc = 0;
    else
        wauc = sum(Wcum(icum,:),2);
    end
    
    % If the 1st or last value is NaN, trim the new X and Y
    if needTrimLast && (isnan(X(end)) || isnan(Y(end)))
        X(end) = [];
        Y(end) = [];
    end
    if needTrimFirst && (isnan(X(1)) || isnan(Y(1)))
        X(1) = [];
        Y(1) = [];
    end
    
    % Have enough data?
    if length(X)<2
        auc = 0;
        return;
    end
    
    % Get area
    auc = 0.5*sum( (X(2:end)-X(1:end-1)).*(Y(2:end)+Y(1:end-1)) );
    auc = abs(auc);
    end
    
    
    function optpt = findoptroc(X,Y,Wcum,scale,cost)
    % Get positive and negative counts
    wP = scale(1)*Wcum(end,1);
    wN = scale(2)*sum(Wcum(end,2:end),2);
    
    % Get the optimal slope
    m = (cost(2,1)-cost(2,2))/(cost(1,2)-cost(1,1)) * wN/wP;
    
    % Find lowest intercept for straight lines drawn through (X,Y) 
    %   using this slope and X axis
    [~,idx] = min(X - Y/m);
    
    % Get the optimal point
    optpt = [X(idx) Y(idx)];
    end
    
    
    function [scores,labels,weights,ncv] = preparedata(scores,labels,weights)
    ncv = 0; % no cross-validation
    if iscell(scores) % Cross-validated scores and labels
        % Scores
        if ~isvector(scores)
            error(message('stats:perfcurve:CVScoresNotVector'));
        end
        scores = scores(:);
        ncv = numel(scores);
        if ncv<2
            error(message('stats:perfcurve:CVScoresTooShort'));
        end
        tfisemp = cellfun(@isempty,scores);
        tfisnum = cellfun(@isnumeric,scores);
        tfisvec = cellfun(@isvector,scores);
        if any(tfisemp) || any(~tfisnum) || any(~tfisvec)
            error(message('stats:perfcurve:CVScoresWithBadElements'));
        end
        nels = cellfun(@numel,scores);
        for icv=1:ncv
            scores{icv} = scores{icv}(:);
        end
        
        % Labels
        if ~iscell(labels)
            error(message('stats:perfcurve:CVLabelsNotCell'));
        end
        if ~isvector(labels)
            error(message('stats:perfcurve:CVLabelsNotVector'));
        end
        if numel(labels)~=ncv
            error(message('stats:perfcurve:CVLabelsWithNonmatchingLength'));
        end
        labels = labels(:);
        ltype = cellfun(@class,labels,'UniformOutput',false);
        for icv=2:ncv
            if ~strcmp(ltype{1},ltype{icv})
                error(message('stats:perfcurve:CVLabelsWithDifferentTypes'));
            end
        end
        nell = zeros(ncv,1);
        for icv=1:ncv
            if strcmp(ltype{icv},'char')
                nell(icv) = size(labels{icv},1);
            else
                labels{icv} = labels{icv}(:);
                nell(icv) = numel(labels{icv});
            end
        end
        
        % Weights
        if isempty(weights)
            weights = cell(ncv,1);
        else
            if ~iscell(weights)
                error(message('stats:perfcurve:CVWeightsNotCell'));
            end
            if ~isvector(weights)
                error(message('stats:perfcurve:CVWeightsNotVector'));
            end
            if numel(weights)~=ncv
                error(message('stats:perfcurve:CVWeightsWithNonmatchingLength'));
            end
            weights = weights(:);
        end
        for icv=1:ncv
            if isempty(weights{icv})
                weights{icv} = ones(nels(icv),1);
            else
                weights{icv} = weights{icv}(:);
            end
        end
        tfisflo = cellfun(@isfloat,weights);
        tfisvec = cellfun(@isvector,weights);
        if any(~tfisflo) || any(~tfisvec)
            error(message('stats:perfcurve:CVWeightsNotNumeric'));
        end
        fneg = @(x) any(x<0);
        fzer = @(x) all(x==0);
        fnan = @(x) all(isnan(x));
        tfisneg = cellfun(fneg,weights);
        tfiszer = cellfun(fzer,weights);
        tfisnan = cellfun(fnan,weights);
        if any(tfisneg) || any(tfiszer) || any(tfisnan)
            error(message('stats:perfcurve:CVWeightsNegativeOrNaN'));
        end
        nelw = cellfun(@numel,weights);
        
        % Check sizes
        if any(nels~=nell) || any(nelw~=nell)
            error(message('stats:perfcurve:CVWeightsNotMatchedToScores'));
        end
        
        % Get rid of observations with zero weights
        for icv=1:ncv
            iszer = weights{icv}==0;
            scores{icv}(iszer) = [];
            labels{icv}(iszer,:) = [];
            weights{icv}(iszer) = [];
        end
    else % Plain scores, labels and weights
        % Scores
        if ~isvector(scores)
            error(message('stats:perfcurve:ScoresNotVector'));
        end
        scores = scores(:);
        ns = numel(scores);
        
        % Labels
        if ~ischar(labels)
            if ~isvector(labels)
                error(message('stats:perfcurve:LabelsNotVector'));
            end
            labels = labels(:);
        end
        nl = size(labels,1);
        if ns~=nl
            error(message('stats:perfcurve:ScoresAndLabelsDoNotMatch'));
        end
        
        % Weights
        if isempty(weights)
            weights = ones(ns,1);
        else
            if ~isfloat(weights)
                error(message('stats:perfcurve:WeightsNotRealNumbers'));
            end
            if ~isvector(weights)
                error(message('stats:perfcurve:WeightsNotVector'));
            end
            if any(weights<0) || all(weights==0)
                error(message('stats:perfcurve:WeightsNegative'));
            end
            if numel(weights)~=ns
                error(message('stats:perfcurve:WeightsNotMatchedToScores'));
            end
            if any(isnan(weights))
                error(message('stats:perfcurve:NaNWeights'));
            end
            weights = weights(:);
            iszer = weights==0;
            scores(iszer) = [];
            labels(iszer,:) = [];
            weights(iszer) = [];
        end
    end
    end
    
    
    % Stack cell arrays of scores used for cross-validation, labels and weights
    % into long vectors preserving type. ncvel is a vector storing the number
    % of elements in every CV piece.
    function [ncvel,scores,labels,weights] = stackcv(scores,labels,weights)
    ncvel = cellfun(@numel,scores);
    scores = vertcat(scores{:});
    weights = vertcat(weights{:});
    if ischar(labels{1})
        labels = char(labels{:});
    else
        labels = vertcat(labels{:});
    end
    end
    
    
    % Compute division indices, and X, Y and T values for either
    % cross-validation or bootstrap samples. scores and Wcum are cell arrays
    % with Nsub elements, one element per subset. wX and wY are weights
    % appropriate for computation of statistics X and Y. For thresholds T,
    % always use the total weight of the subset for any threshold.
    function [X,wX,Y,wY,T,wT] = xyt(xVals,tVals,scores,Wcum,afx,afwx,afy,afwy)
    % Init
    if     ~isempty(xVals)
        M = length(xVals);
    else
    %elseif ~isempty(tVals)
        M = length(tVals);
    end
    nsub = length(scores);
    X = zeros(M,nsub);
    Y = zeros(M,nsub);
    wX = zeros(M,nsub);
    wY = zeros(M,nsub);
    
    % Loop through subsets
    T = [];
    wT = [];
    computeT = nargout>4 && ~isempty(xVals);
    if computeT
        T = zeros(M,nsub);
        wT = zeros(M,nsub);
        for isub=1:nsub
            [X(:,isub),wX(:,isub),Y(:,isub),wY(:,isub),T(:,isub),wT(:,isub)] = ...
                xytOneSample(xVals,tVals,scores{isub},Wcum{isub},afx,afwx,afy,afwy);
        end
    else
        for isub=1:nsub
            [X(:,isub),wX(:,isub),Y(:,isub),wY(:,isub)] = ...
                xytOneSample(xVals,tVals,scores{isub},Wcum{isub},afx,afwx,afy,afwy);
        end
    end
    end
    
    
    % Returns one sample of XYT values with corresponding weights.
    % Weights for T averaging are set to the cumulative weight of the sample.
    function [X,wX,Y,wY,T,wT] = xytOneSample(xVals,tVals,scores,Wcum,afx,afwx,afy,afwy)
    % Notify the user once that the class is missing in one of the folds.
    if ~all(sum(Wcum,1)) && notifyClassAbsent()
        warning(message('stats:perfcurve:SubSampleWithMissingClasses'));
    end
    
    % Get divisions
    if     ~isempty(tVals) % perform threshold averaging
        div = tdiv(tVals,scores,false);
    else
    %elseif ~isempty(xVals) % perform vertical averaging
        div = xdiv(xVals,afx,Wcum,false);
    end
    
    % Get X and Y
    [X,wX,tpX,fpX] = Xvalues(div,Wcum,afx,afwx);
    [Y,wY] = Yvalues(tpX,fpX,Wcum,afy,afwy);
    
    % Get T
    T = [];
    wT = [];
    computeT = nargout>4 && ~isempty(xVals);
    if computeT
        if     ~isempty(xVals)
            M = length(xVals);
        elseif ~isempty(tVals)
            M = length(tVals);
        end
        T = thresholds(div,scores,false);
        wT = repmat(sum(Wcum(end,:),2),M,1);
    end
    end
    
    
    % Returns CV confidence interval for every row of Xcv.
    function ci = cvci(X,Xcv,wXcv,alpha)
    ci = [];
    ncv = size(Xcv,2);
    if ncv>0
        e = wstd(X,Xcv,wXcv)/sqrt(ncv);
        ea = e*norminv(1-0.5*alpha,0,1);
        ci = [X-ea X+ea];
    end
    end
    
    
    % Returns one bootstrap replica specified by indices idx.
    % Weights are not needed here because they are used to generate replicas
    % and must be fed into bootci directly.
    function out = oneBootXYT(idx,scores,C,xVals,tVals,processNaN,afx,afy)
    idx = sort(idx);
    [Ccum,subscores] = makeccum(C(idx,:),scores(idx),processNaN);
    if     isempty(xVals) % TA
        [X,~,Y] = xytOneSample(xVals,tVals,subscores,Ccum,afx,[],afy,[]);
        out = [X Y];
    else
    %elseif isempty(tVals) % VA
        [~,~,Y,~,T,~] = xytOneSample(xVals,tVals,subscores,Ccum,afx,[],afy,[]);
        out = [T Y];
    end
    end
    
    
    function out = oneBootAUC(idx,scores,C,XorTrange,mode,processNaN,afx,afy)
    idx = sort(idx);
    [Ccum,subscores] = makeccum(C(idx,:),scores(idx),processNaN);
    out = AUC(XorTrange,mode,subscores,Ccum,afx,afy);
    end
    
    
    % Compute weighted standard deviation given mean M, 
    % measurements X and their weights W.
    % X and W are NxK matrices for N measurements and K folds.
    function E = wstd(M,X,W)
    tfnan = isnan(X);
    X(tfnan) = 0;
    W(tfnan) = 0;
    tfzeroW = ~any(W,2);
    if any(tfzeroW)
        warning(message('stats:perfcurve:CVfoldsWithZeroWeights'));
    end
    E = zeros(length(M),1);
    i = ~tfzeroW;
    E(i) = sqrt(sum(W(i,:).*bsxfun(@minus,X(i,:),M(i)).^2,2) ./ sum(W(i,:),2));
    end
    
    
    % Notify the user if a class is absent?
    function tf = notifyClassAbsent(reset)
    persistent notified;
    if nargin>0 && reset
        tf = true;
        notified = false;
    else
        tf = ~notified;
        notified = true;
    end
    end
    

      

  • 相关阅读:
    HOJ 2139 Spiderman's workout(动态规划)
    FZU 2107 Hua Rong Dao(dfs)
    Java 第十一届 蓝桥杯 省模拟赛 计算机存储中有多少字节
    Java 第十一届 蓝桥杯 省模拟赛 计算机存储中有多少字节
    Java 第十一届 蓝桥杯 省模拟赛 计算机存储中有多少字节
    Java 第十一届 蓝桥杯 省模拟赛 合法括号序列
    Java 第十一届 蓝桥杯 省模拟赛 合法括号序列
    Java 第十一届 蓝桥杯 省模拟赛 合法括号序列
    Java 第十一届 蓝桥杯 省模拟赛 无向连通图最少包含多少条边
    Java 第十一届 蓝桥杯 省模拟赛 无向连通图最少包含多少条边
  • 原文地址:https://www.cnblogs.com/sword-/p/9357224.html
Copyright © 2011-2022 走看看