Daniel@0: function Perf = mlr_test(W, test_k, Xtrain, Ytrain, Xtest, Ytest, Testnorm)
Daniel@0: %   Perf = mlr_test(W, test_k, Xtrain, Ytrain, Xtest, Ytest)
Daniel@0: %
Daniel@0: %       W       = d-by-d positive semi-definite matrix
Daniel@0: %       test_k  = vector of k-values to use for KNN/Prec@k/NDCG
Daniel@0: %       Xtrain  = d-by-n matrix of training data
Daniel@0: %       Ytrain  = n-by-1 vector of training labels
Daniel@0: %                   OR
Daniel@0: %                 n-by-2 cell array where
Daniel@0: %                   Y{q,1} contains relevant indices (in 1..n) for point q
Daniel@0: %                   Y{q,2} contains irrelevant indices (in 1..n) for point q
Daniel@0: %       Xtest   = d-by-m matrix of testing data
Daniel@0: %       Ytest   = m-by-1 vector of training labels, or m-by-2 cell array 
Daniel@0: %               
Daniel@0: %
Daniel@0: %   The output structure Perf contains the mean score for:
Daniel@0: %       AUC, KNN, Prec@k, MAP, MRR, NDCG,
Daniel@0: %   as well as the effective dimensionality of W, and
Daniel@0: %   the best-performing k-value for KNN, Prec@k, and NDCG.
Daniel@0: %
Daniel@0: 
Daniel@0: %     addpath('cuttingPlane', 'distance', 'feasible', 'initialize', 'loss', ...
Daniel@0: %             'metricPsi', 'regularize', 'separationOracle', 'util');
Daniel@0: 
Daniel@0:     Perf        = struct(                       ...
Daniel@0:                             'AUC',      [],     ...
Daniel@0:                             'KNN',      [],     ...
Daniel@0:                             'PrecAtK',  [],     ...
Daniel@0:                             'MAP',      [],     ...
Daniel@0:                             'MRR',      [],     ...
Daniel@0:                             'NDCG',     [],     ...
Daniel@0:                             'dimensionality',   [],     ...
Daniel@0:                             'KNNk',     [],     ...
Daniel@0:                             'PrecAtKk', [],     ...
Daniel@0:                             'NDCGk',    []     ...
Daniel@0:                 );
Daniel@0: 
Daniel@0:     [d, nTrain, nKernel] = size(Xtrain);
Daniel@0:     nTest       = length(Ytest);
Daniel@0:     test_k      = min(test_k, nTrain);
Daniel@0: 
Daniel@0:     if nargin < 7
Daniel@0:         Testnorm = [];
Daniel@0:     end
Daniel@0: 
Daniel@0:     % Compute dimensionality of the learned metric
Daniel@0:     Perf.dimensionality = mlr_test_dimension(W, nTrain, nKernel);
Daniel@0: 
Daniel@0:     
Daniel@0:     % Build the distance matrix
Daniel@0:     [D, I] = mlr_test_distance(W, Xtrain, Xtest, Testnorm);
Daniel@0: 
Daniel@0: 
Daniel@0:     % Compute label agreement
Daniel@0:     if ~iscell(Ytest)
Daniel@0:         Labels  = Ytrain(I);
Daniel@0:         Agree   = bsxfun(@eq, Ytest', Labels); 
Daniel@0: 
Daniel@0:         % We only compute KNN error if Y are labels
Daniel@0:         [Perf.KNN, Perf.KNNk] = mlr_test_knn(Labels, Ytest, test_k);
Daniel@0:     else
Daniel@0:         Agree   = zeros(nTrain, nTest);
Daniel@0:         for i = 1:nTest
Daniel@0:             Agree(:,i) = ismember(I(:,i), Ytest{i,1});
Daniel@0:         end
Daniel@0:         Agree = reduceAgreement(Agree);
Daniel@0:     end
Daniel@0: 
Daniel@0:     % Compute AUC score
Daniel@0:     Perf.AUC    = mlr_test_auc(Agree);
Daniel@0: 
Daniel@0:     % Compute MAP score
Daniel@0:     Perf.MAP    = mlr_test_map(Agree);
Daniel@0: 
Daniel@0:     % Compute MRR score
Daniel@0:     Perf.MRR    = mlr_test_mrr(Agree);
Daniel@0: 
Daniel@0:     % Compute prec@k
Daniel@0:     [Perf.PrecAtK, Perf.PrecAtKk] = mlr_test_preck(Agree, test_k);
Daniel@0: 
Daniel@0:     % Compute NDCG score
Daniel@0:     [Perf.NDCG, Perf.NDCGk] = mlr_test_ndcg(Agree, test_k);
Daniel@0: 
Daniel@0: end
Daniel@0: 
Daniel@0: 
Daniel@0: function [D,I] = mlr_test_distance(W, Xtrain, Xtest, Testnorm)
Daniel@0: 
Daniel@0:     % CASES:
Daniel@0:     %   Raw:                        W = []
Daniel@0:     
Daniel@0:     %   Linear, full:               W = d-by-d
Daniel@0:     %   Single Kernel, full:        W = n-by-n
Daniel@0:     %   MKL, full:                  W = n-by-n-by-m
Daniel@0: 
Daniel@0:     %   Linear, diagonal:           W = d-by-1
Daniel@0:     %   Single Kernel, diagonal:    W = n-by-1
Daniel@0:     %   MKL, diag:                  W = n-by-m
Daniel@0:     %   MKL, diag-off-diag:         W = m-by-m-by-n
Daniel@0:     
Daniel@0:     [d, nTrain, nKernel] = size(Xtrain);
Daniel@0:     nTest = size(Xtest, 2);
Daniel@0: 
Daniel@0:     if isempty(W)
Daniel@0:         % W = []  => native euclidean distances
Daniel@0:         D = mlr_test_distance_raw(Xtrain, Xtest, Testnorm);
Daniel@0: 
Daniel@0:     elseif size(W,1) == d && size(W,2) == d
Daniel@0:         % We're in a full-projection case
Daniel@0:         D = setDistanceFullMKL([Xtrain Xtest], W, nTrain + (1:nTest), 1:nTrain);
Daniel@0: 
Daniel@0:     elseif size(W,1) == d && size(W,2) == nKernel
Daniel@0:         % We're in a simple diagonal case
Daniel@0:         D = setDistanceDiagMKL([Xtrain Xtest], W, nTrain + (1:nTest), 1:nTrain);
Daniel@0: 
Daniel@0:     elseif size(W,1) == nKernel && size(W,2) == nKernel && size(W,3) == nTrain
Daniel@0:         % We're in DOD mode
Daniel@0:         D = setDistanceDODMKL([Xtrain Xtest], W, nTrain + (1:nTest), 1:nTrain);
Daniel@0: 
Daniel@0:     else
Daniel@0:         % Error?
Daniel@0:         error('Cannot determine metric mode.');
Daniel@0: 
Daniel@0:     end
Daniel@0:     
Daniel@0:     D       = full(D(1:nTrain, nTrain + (1:nTest)));
Daniel@0:     [v,I]   = sort(D, 1);
Daniel@0: end
Daniel@0: 
Daniel@0: 
Daniel@0: 
Daniel@0: function dimension = mlr_test_dimension(W, nTrain, nKernel)
Daniel@0: 
Daniel@0:     % CASES:
Daniel@0:     %   Raw:                        W = []
Daniel@0:     
Daniel@0:     %   Linear, full:               W = d-by-d
Daniel@0:     %   Single Kernel, full:        W = n-by-n
Daniel@0:     %   MKL, full:                  W = n-by-n-by-m
Daniel@0: 
Daniel@0:     %   Linear, diagonal:           W = d-by-1
Daniel@0:     %   Single Kernel, diagonal:    W = n-by-1
Daniel@0:     %   MKL, diag:                  W = n-by-m
Daniel@0:     %   MKL, diag-off-diag:         W = m-by-m-by-n
Daniel@0:     
Daniel@0:     
Daniel@0:     if size(W,2) == nTrain
Daniel@0:         dim = [];
Daniel@0:         for i = 1:nKernel
Daniel@0:             [v,d]   = eig(0.5 * (W(:,:,i) + W(:,:,i)'));
Daniel@0:             dim     = [dim ; abs(real(diag(d)))];
Daniel@0:         end
Daniel@0:     else
Daniel@0:         dim       = W(:);
Daniel@0:     end
Daniel@0: 
Daniel@0:     cd      = cumsum(dim) / sum(dim);
Daniel@0:     dimension = find(cd >= 0.95, 1);
Daniel@0:     if isempty(dimension)
Daniel@0:         dimension = 0;
Daniel@0:     end
Daniel@0: end
Daniel@0: 
Daniel@0: function [NDCG, NDCGk] = mlr_test_ndcg(Agree, test_k)
Daniel@0: 
Daniel@0:     nTrain = size(Agree, 1);
Daniel@0: 
Daniel@0:     Discount        = zeros(1, nTrain);
Daniel@0:     Discount(1:2)   = 1;
Daniel@0: 
Daniel@0:     NDCG   = -Inf;
Daniel@0:     NDCGk  = 0;
Daniel@0:     for k = test_k
Daniel@0:         
Daniel@0:         Discount(3:k)   = 1 ./ log2(3:k);
Daniel@0:         Discount        = Discount / sum(Discount);
Daniel@0: 
Daniel@0:         b = mean(Discount * Agree);
Daniel@0:         if b > NDCG
Daniel@0:             NDCG = b;
Daniel@0:             NDCGk = k;
Daniel@0:         end
Daniel@0:     end
Daniel@0: end
Daniel@0: 
Daniel@0: function [PrecAtK, PrecAtKk] = mlr_test_preck(Agree, test_k)
Daniel@0: 
Daniel@0:     PrecAtK        = -Inf;
Daniel@0:     PrecAtKk       = 0;
Daniel@0:     for k = test_k
Daniel@0:         b   = mean( mean( Agree(1:k, :), 1 ) );
Daniel@0:         if b > PrecAtK
Daniel@0:             PrecAtK = b;
Daniel@0:             PrecAtKk = k;
Daniel@0:         end
Daniel@0:     end
Daniel@0: end
Daniel@0: 
Daniel@0: function [KNN, KNNk] = mlr_test_knn(Labels, Ytest, test_k)
Daniel@0: 
Daniel@0:     KNN        = -Inf;
Daniel@0:     KNNk       = 0;
Daniel@0:     for k = test_k
Daniel@0:         b   = mean( mode( Labels(1:k,:), 1 ) == Ytest');
Daniel@0:         if b > KNN
Daniel@0:             KNN    = b;
Daniel@0:             KNNk   = k;
Daniel@0:         end
Daniel@0:     end
Daniel@0: end
Daniel@0: 
Daniel@0: function MAP = mlr_test_map(Agree);
Daniel@0: 
Daniel@0:     nTrain      = size(Agree, 1);
Daniel@0:     MAP         = bsxfun(@ldivide, (1:nTrain)', cumsum(Agree, 1));
Daniel@0:     MAP         = mean(sum(MAP .* Agree, 1)./ sum(Agree, 1));
Daniel@0: end
Daniel@0: 
Daniel@0: function MRR = mlr_test_mrr(Agree);
Daniel@0: 
Daniel@0:         nTest = size(Agree, 2);
Daniel@0:         MRR        = 0;
Daniel@0:         for i = 1:nTest
Daniel@0:             MRR    = MRR  + (1 / find(Agree(:,i), 1));
Daniel@0:         end
Daniel@0:         MRR        = MRR / nTest;
Daniel@0: end
Daniel@0: 
Daniel@0: function AUC = mlr_test_auc(Agree)
Daniel@0: 
Daniel@0:     TPR             = cumsum(Agree,     1);
Daniel@0:     FPR             = cumsum(~Agree,    1);
Daniel@0: 
Daniel@0:     numPos          = TPR(end,:);
Daniel@0:     numNeg          = FPR(end,:);
Daniel@0: 
Daniel@0:     TPR             = mean(bsxfun(@rdivide, TPR, numPos),2);
Daniel@0:     FPR             = mean(bsxfun(@rdivide, FPR, numNeg),2);
Daniel@0:     AUC             = diff([0 FPR']) * TPR;
Daniel@0: end
Daniel@0: 
Daniel@0: 
Daniel@0: function D = mlr_test_distance_raw(Xtrain, Xtest, Testnorm)
Daniel@0: 
Daniel@0:     [d, nTrain, nKernel] = size(Xtrain);
Daniel@0:     nTest = size(Xtest, 2);
Daniel@0: 
Daniel@0:     if isempty(Testnorm)
Daniel@0:         % Not in kernel mode, compute distances directly
Daniel@0:         D = 0;
Daniel@0:         for i = 1:nKernel
Daniel@0:             D = D + setDistanceDiag([Xtrain(:,:,i) Xtest(:,:,i)], ones(d,1), ...
Daniel@0:                                     nTrain + (1:nTest), 1:nTrain);
Daniel@0:         end
Daniel@0:     else
Daniel@0:         % We are in kernel mode
Daniel@0:         D = sparse(nTrain + nTest, nTrain + nTest);
Daniel@0:         for i = 1:nKernel
Daniel@0:             Trainnorm = diag(Xtrain(:,:,i));
Daniel@0:             D(1:nTrain, nTrain + (1:nTest)) = D(1:nTrain, nTrain + (1:nTest)) ...
Daniel@0:                 +  bsxfun(@plus, Trainnorm, bsxfun(@plus, Testnorm(:,i)', -2 * Xtest(:,:,i)));
Daniel@0:         end
Daniel@0:     end
Daniel@0: end
Daniel@0: 
Daniel@0: function A = reduceAgreement(Agree)
Daniel@0:     nPos = sum(Agree,1);
Daniel@0:     nNeg = sum(~Agree,1);
Daniel@0: 
Daniel@0:     goodI = find(nPos > 0 & nNeg > 0);
Daniel@0:     A = Agree(:,goodI);
Daniel@0: end