Daniel@0: function [Y, Loss] = separationOracleNDCG(q, D, pos, neg, k)
Daniel@0: %
Daniel@0: %   [Y,Loss]  = separationOracleNDCG(q, D, pos, neg, k)
Daniel@0: %
Daniel@0: %   q   = index of the query point
Daniel@0: %   D   = the current distance matrix
Daniel@0: %   pos = indices of relevant results for q
Daniel@0: %   neg = indices of irrelevant results for q
Daniel@0: %   k   = length of the list to consider 
Daniel@0: %
Daniel@0: %   Y is a permutation 1:n corresponding to the maximally
Daniel@0: %   violated constraint
Daniel@0: %
Daniel@0: %   Loss is the loss for Y, in this case, 1-NDCG(Y)
Daniel@0: 
Daniel@0: 
Daniel@0:     % First, sort the documents in descending order of W'Phi(q,x)
Daniel@0:     % Phi = - (X(q) - X(x)) * (X(q) - X(x))'
Daniel@0:     
Daniel@0:     % Sort the positive documents
Daniel@0:     ScorePos        = - D(pos, q);
Daniel@0:     [Vpos, Ipos]    = sort(full(ScorePos'), 'descend');
Daniel@0:     Ipos            = pos(Ipos);
Daniel@0:     
Daniel@0:     % Sort the negative documents
Daniel@0:     ScoreNeg        = - D(neg, q);
Daniel@0:     [Vneg, Ineg]    = sort(full(ScoreNeg'), 'descend');
Daniel@0:     Ineg            = neg(Ineg);
Daniel@0: 
Daniel@0:     % Now, solve the DP for the interleaving
Daniel@0: 
Daniel@0:     numPos  = length(pos);
Daniel@0:     numNeg  = length(neg);
Daniel@0:     n       = numPos + numNeg;
Daniel@0: 
Daniel@0:     % From Chakrabarti (KDD08)
Daniel@0:     k       = min(k, numPos);
Daniel@0: 
Daniel@0:     cVneg   = cumsum(Vneg);
Daniel@0: 
Daniel@0:     Discount = zeros(k, 1);
Daniel@0:     Discount(1:2) = 1;
Daniel@0:     Discount(3:k) = 1./ log2(3:k);
Daniel@0: 
Daniel@0:     DCGstar         = sum(Discount);
Daniel@0: 
Daniel@0: 
Daniel@0:     % Pre-compute the loss table
Daniel@0:     LossTab     = padarray(  hankel(- Discount / DCGstar), ...
Daniel@0:                             max(0, [numNeg numPos] - k), 0, 'post');
Daniel@0:     if sum(size(LossTab) > [numNeg, numPos])
Daniel@0:         LossTab     = LossTab(1:numNeg, 1:numPos);
Daniel@0:     end
Daniel@0: 
Daniel@0:     % 2010-01-17 09:13:41 by Brian McFee <bmcfee@cs.ucsd.edu>
Daniel@0:     % initialize the score table
Daniel@0: 
Daniel@0:     pcVneg = [0 cVneg];
Daniel@0:     % Pre-compute cellScore
Daniel@0:     cellValue = bsxfun(@times, Vpos / (numPos * numNeg), numNeg - 2 * ((1:numNeg)-1)');
Daniel@0:     cellValue = bsxfun(@plus, (2 * pcVneg(1:numNeg) - cVneg(end))' / (numPos * numNeg), cellValue);
Daniel@0:     cellValue = cellValue + LossTab;
Daniel@0: 
Daniel@0:     S       = zeros(numNeg, numPos);
Daniel@0:     P       = zeros(numNeg, numPos);
Daniel@0:     
Daniel@0:     % Initialize first column
Daniel@0:     P(:,1) = 1;
Daniel@0:     S(:,1) = cellValue(:,1);
Daniel@0: 
Daniel@0:     % Initialize first row
Daniel@0:     P(1,:) = 1;
Daniel@0:     S(1,:) = cumsum(cellValue(1,:));
Daniel@0:     
Daniel@0:     % For the rest, use the recurrence
Daniel@0: 
Daniel@0:     for g = 2:numPos
Daniel@0:         [m, pointer]    = cummax(S(:,g-1));
Daniel@0:         P(:,g)          = pointer;
Daniel@0:         S(:,g)          = m' + cellValue(:,g);
Daniel@0:     end
Daniel@0: 
Daniel@0:     % Now reconstruct the permutation from the DP table
Daniel@0:     Y           = nan * ones(n,1);
Daniel@0:     [m,p]       = max(S(:,numPos));
Daniel@0:     
Daniel@0:     Loss        = 1 + LossTab(p,numPos);
Daniel@0:     
Daniel@0:     NegsBefore      = zeros(numPos,1);
Daniel@0:     NegsBefore(numPos)  = p-1;
Daniel@0: 
Daniel@0:     for a = numPos:-1:2
Daniel@0:         p               = P(p,a);
Daniel@0:         NegsBefore(a-1) = p-1;
Daniel@0:         Loss            = Loss + LossTab(p,a-1);
Daniel@0:     end
Daniel@0:     Y((1:numPos)' + NegsBefore)     = Ipos;
Daniel@0:     Y(isnan(Y))                     = Ineg;
Daniel@0: 
Daniel@0: end