wolffd@0: function Perf = mlr_test(W, test_k, Xtrain, Ytrain, Xtest, Ytest)
wolffd@0: %   Perf = mlr_test(W, test_k, Xtrain, Ytrain, Xtest, Ytest)
wolffd@0: %
wolffd@0: %       W       = d-by-d positive semi-definite matrix
wolffd@0: %       test_k  = vector of k-values to use for KNN/Prec@k/NDCG
wolffd@0: %       Xtrain  = d-by-n matrix of training data
wolffd@0: %       Ytrain  = n-by-1 vector of training labels
wolffd@0: %                   OR
wolffd@0: %                 n-by-2 cell array where
wolffd@0: %                   Y{q,1} contains relevant indices (in 1..n) for point q
wolffd@0: %                   Y{q,2} contains irrelevant indices (in 1..n) for point q
wolffd@0: %       Xtest   = d-by-m matrix of testing data
wolffd@0: %       Ytest   = m-by-1 vector of training labels, or m-by-2 cell array 
wolffd@0: %               
wolffd@0: %
wolffd@0: %   The output structure Perf contains the mean score for:
wolffd@0: %       AUC, KNN, Prec@k, MAP, MRR, NDCG,
wolffd@0: %   as well as the effective dimensionality of W, and
wolffd@0: %   the best-performing k-value for KNN, Prec@k, and NDCG.
wolffd@0: %
wolffd@0: 
wolffd@0:     Perf        = struct(                       ...
wolffd@0:                             'AUC',      [],     ...
wolffd@0:                             'KNN',      [],     ...
wolffd@0:                             'PrecAtK',  [],     ...
wolffd@0:                             'MAP',      [],     ...
wolffd@0:                             'MRR',      [],     ...
wolffd@0:                             'NDCG',     [],     ...
wolffd@0:                             'dimensionality',   [],     ...
wolffd@0:                             'KNNk',     [],     ...
wolffd@0:                             'PrecAtKk', [],     ...
wolffd@0:                             'NDCGk',    []     ...
wolffd@0:                 );
wolffd@0: 
wolffd@0:     [d, nTrain, nKernel] = size(Xtrain);
wolffd@0:     % Compute dimensionality of the learned metric
wolffd@0:     Perf.dimensionality = mlr_test_dimension(W, nTrain, nKernel);
wolffd@0:     test_k      = min(test_k, nTrain);
wolffd@0: 
wolffd@0:     if nargin > 5
wolffd@0:         % Knock out the points with no labels
wolffd@0:         if ~iscell(Ytest)
wolffd@0:             Ibad                = find(isnan(Ytrain));
wolffd@0:             Xtrain(:,Ibad,:)    = inf;
wolffd@0:         end
wolffd@0: 
wolffd@0:         % Build the distance matrix
wolffd@0:         [D, I] = mlr_test_distance(W, Xtrain, Xtest);
wolffd@0:     else
wolffd@0:         % Leave-one-out validation
wolffd@0: 
wolffd@0:         if nargin > 4 
wolffd@0:             % In this case, Xtest is a subset of training indices to test on
wolffd@0:             testRange = Xtest;
wolffd@0:         else
wolffd@0:             testRange = 1:nTrain;
wolffd@0:         end
wolffd@0:         Xtest       = Xtrain(:,testRange,:);
wolffd@0:         Ytest       = Ytrain(testRange);
wolffd@0: 
wolffd@0:         % compute self-distance
wolffd@0:         [D, I]  = mlr_test_distance(W, Xtrain, Xtest);
wolffd@0:         % clear out the self-link (distance = 0)
wolffd@0:         I       = I(2:end,:);
wolffd@0:         D       = D(2:end,:);
wolffd@0:     end
wolffd@0:     
wolffd@0:     nTest       = length(Ytest);
wolffd@0: 
wolffd@0:     % Compute label agreement
wolffd@0:     if ~iscell(Ytest)
wolffd@0:         % First, knock out the points with no label
wolffd@0:         Labels  = Ytrain(I);
wolffd@0:         Agree   = bsxfun(@eq, Ytest', Labels); 
wolffd@0: 
wolffd@0:         % We only compute KNN error if Y are labels
wolffd@0:         [Perf.KNN, Perf.KNNk] = mlr_test_knn(Labels, Ytest, test_k);
wolffd@0:     else
wolffd@0:         if nargin > 5
wolffd@0:             Agree   = zeros(nTrain, nTest);
wolffd@0:         else
wolffd@0:             Agree   = zeros(nTrain-1, nTest);
wolffd@0:         end
wolffd@0:         for i = 1:nTest
wolffd@0:             Agree(:,i) = ismember(I(:,i), Ytest{i,1});
wolffd@0:         end
wolffd@0: 
wolffd@0:         Agree = reduceAgreement(Agree);
wolffd@0:     end
wolffd@0: 
wolffd@0:     % Compute AUC score
wolffd@0:     Perf.AUC    = mlr_test_auc(Agree);
wolffd@0: 
wolffd@0:     % Compute MAP score
wolffd@0:     Perf.MAP    = mlr_test_map(Agree);
wolffd@0: 
wolffd@0:     % Compute MRR score
wolffd@0:     Perf.MRR    = mlr_test_mrr(Agree);
wolffd@0: 
wolffd@0:     % Compute prec@k
wolffd@0:     [Perf.PrecAtK, Perf.PrecAtKk] = mlr_test_preck(Agree, test_k);
wolffd@0: 
wolffd@0:     % Compute NDCG score
wolffd@0:     [Perf.NDCG, Perf.NDCGk] = mlr_test_ndcg(Agree, test_k);
wolffd@0: 
wolffd@0: end
wolffd@0: 
wolffd@0: 
wolffd@0: function [D,I] = mlr_test_distance(W, Xtrain, Xtest)
wolffd@0: 
wolffd@0:     % CASES:
wolffd@0:     %   Raw:                        W = []
wolffd@0:     
wolffd@0:     %   Linear, full:               W = d-by-d
wolffd@0:     %   Single Kernel, full:        W = n-by-n
wolffd@0:     %   MKL, full:                  W = n-by-n-by-m
wolffd@0: 
wolffd@0:     %   Linear, diagonal:           W = d-by-1
wolffd@0:     %   Single Kernel, diagonal:    W = n-by-1
wolffd@0:     %   MKL, diag:                  W = n-by-m
wolffd@0:     %   MKL, diag-off-diag:         W = m-by-m-by-n
wolffd@0:     
wolffd@0:     [d, nTrain, nKernel] = size(Xtrain);
wolffd@0:     nTest = size(Xtest, 2);
wolffd@0: 
wolffd@0:     if isempty(W)
wolffd@0:         % W = []  => native euclidean distances
wolffd@0:         D = mlr_test_distance_raw(Xtrain, Xtest);
wolffd@0: 
wolffd@0:     elseif size(W,1) == d && size(W,2) == d
wolffd@0:         % We're in a full-projection case
wolffd@0:         D = setDistanceFullMKL([Xtrain Xtest], W, nTrain + (1:nTest), 1:nTrain);
wolffd@0: 
wolffd@0:     elseif size(W,1) == d && size(W,2) == nKernel
wolffd@0:         % We're in a simple diagonal case
wolffd@0:         D = setDistanceDiagMKL([Xtrain Xtest], W, nTrain + (1:nTest), 1:nTrain);
wolffd@0: 
wolffd@0:     else
wolffd@0:         % Error?
wolffd@0:         error('Cannot determine metric mode.');
wolffd@0: 
wolffd@0:     end
wolffd@0:     
wolffd@0:     D       = full(D(1:nTrain, nTrain + (1:nTest)));
wolffd@0:     [v,I]   = sort(D, 1);
wolffd@0: end
wolffd@0: 
wolffd@0: 
wolffd@0: 
wolffd@0: function dimension = mlr_test_dimension(W, nTrain, nKernel)
wolffd@0: 
wolffd@0:     % CASES:
wolffd@0:     %   Raw:                        W = []
wolffd@0:     
wolffd@0:     %   Linear, full:               W = d-by-d
wolffd@0:     %   Single Kernel, full:        W = n-by-n
wolffd@0:     %   MKL, full:                  W = n-by-n-by-m
wolffd@0: 
wolffd@0:     %   Linear, diagonal:           W = d-by-1
wolffd@0:     %   Single Kernel, diagonal:    W = n-by-1
wolffd@0:     %   MKL, diag:                  W = n-by-m
wolffd@0:     %   MKL, diag-off-diag:         W = m-by-m-by-n
wolffd@0:     
wolffd@0:     
wolffd@0:     if size(W,1) == size(W,2)
wolffd@0:         dim = [];
wolffd@0:         for i = 1:nKernel
wolffd@0:             [v,d]   = eig(0.5 * (W(:,:,i) + W(:,:,i)'));
wolffd@0:             dim     = [dim ; abs(real(diag(d)))];
wolffd@0:         end
wolffd@0:     else
wolffd@0:         dim       = W(:);
wolffd@0:     end
wolffd@0: 
wolffd@0:     cd      = cumsum(dim) / sum(dim);
wolffd@0:     dimension = find(cd >= 0.95, 1);
wolffd@0:     if isempty(dimension)
wolffd@0:         dimension = 0;
wolffd@0:     end
wolffd@0: end
wolffd@0: 
wolffd@0: function [NDCG, NDCGk] = mlr_test_ndcg(Agree, test_k)
wolffd@0: 
wolffd@0:     nTrain = size(Agree, 1);
wolffd@0: 
wolffd@0:     Discount        = zeros(1, nTrain);
wolffd@0:     Discount(1:2)   = 1;
wolffd@0: 
wolffd@0:     NDCG   = -Inf;
wolffd@0:     NDCGk  = 0;
wolffd@0:     for k = test_k
wolffd@0:         
wolffd@0:         Discount(3:k)   = 1 ./ log2(3:k);
wolffd@0:         Discount        = Discount / sum(Discount);
wolffd@0: 
wolffd@0:         b = mean(Discount * Agree);
wolffd@0:         if b > NDCG
wolffd@0:             NDCG = b;
wolffd@0:             NDCGk = k;
wolffd@0:         end
wolffd@0:     end
wolffd@0: end
wolffd@0: 
wolffd@0: function [PrecAtK, PrecAtKk] = mlr_test_preck(Agree, test_k)
wolffd@0: 
wolffd@0:     PrecAtK        = -Inf;
wolffd@0:     PrecAtKk       = 0;
wolffd@0:     for k = test_k
wolffd@0:         b   = mean( mean( Agree(1:k, :), 1 ) );
wolffd@0:         if b > PrecAtK
wolffd@0:             PrecAtK = b;
wolffd@0:             PrecAtKk = k;
wolffd@0:         end
wolffd@0:     end
wolffd@0: end
wolffd@0: 
wolffd@0: function [KNN, KNNk] = mlr_test_knn(Labels, Ytest, test_k)
wolffd@0: 
wolffd@0:     KNN        = -Inf;
wolffd@0:     KNNk       = 0;
wolffd@0:     for k = test_k
wolffd@0:         % FIXME:  2012-02-07 16:51:59 by Brian McFee <bmcfee@cs.ucsd.edu>
wolffd@0:         %   fix these to discount nans 
wolffd@0: 
wolffd@0:         b   = mean( mode( Labels(1:k,:), 1 ) == Ytest');
wolffd@0:         if b > KNN
wolffd@0:             KNN    = b;
wolffd@0:             KNNk   = k;
wolffd@0:         end
wolffd@0:     end
wolffd@0: end
wolffd@0: 
wolffd@0: function MAP = mlr_test_map(Agree);
wolffd@0: 
wolffd@0:     nTrain      = size(Agree, 1);
wolffd@0:     MAP         = bsxfun(@ldivide, (1:nTrain)', cumsum(Agree, 1));
wolffd@0:     MAP         = mean(sum(MAP .* Agree, 1)./ sum(Agree, 1));
wolffd@0: end
wolffd@0: 
wolffd@0: function MRR = mlr_test_mrr(Agree);
wolffd@0: 
wolffd@0:         nTest = size(Agree, 2);
wolffd@0:         MRR        = 0;
wolffd@0:         for i = 1:nTest
wolffd@0:             MRR    = MRR  + (1 / find(Agree(:,i), 1));
wolffd@0:         end
wolffd@0:         MRR        = MRR / nTest;
wolffd@0: end
wolffd@0: 
wolffd@0: function AUC = mlr_test_auc(Agree)
wolffd@0: 
wolffd@0:     TPR             = cumsum(Agree,     1);
wolffd@0:     FPR             = cumsum(~Agree,    1);
wolffd@0: 
wolffd@0:     numPos          = TPR(end,:);
wolffd@0:     numNeg          = FPR(end,:);
wolffd@0: 
wolffd@0:     TPR             = mean(bsxfun(@rdivide, TPR, numPos),2);
wolffd@0:     FPR             = mean(bsxfun(@rdivide, FPR, numNeg),2);
wolffd@0:     AUC             = diff([0 FPR']) * TPR;
wolffd@0: end
wolffd@0: 
wolffd@0: 
wolffd@0: function D = mlr_test_distance_raw(Xtrain, Xtest)
wolffd@0: 
wolffd@0:     [d, nTrain, nKernel] = size(Xtrain);
wolffd@0:     nTest = size(Xtest, 2);
wolffd@0: 
wolffd@0:         % Not in kernel mode, compute distances directly
wolffd@0:         D = 0;
wolffd@0:         for i = 1:nKernel
wolffd@0:             D = D + setDistanceDiag([Xtrain(:,:,i) Xtest(:,:,i)], ones(d,1), ...
wolffd@0:                                     nTrain + (1:nTest), 1:nTrain);
wolffd@0:         end
wolffd@0: end
wolffd@0: 
wolffd@0: function A = reduceAgreement(Agree)
wolffd@0:     nPos = sum(Agree,1);
wolffd@0:     nNeg = sum(~Agree,1);
wolffd@0: 
wolffd@0:     goodI = find(nPos > 0 & nNeg > 0);
wolffd@0:     A = Agree(:,goodI);
wolffd@0: end