function optA = linear_metric_learning_cvx(data, S, D)
% function optA = linear_metric_learning_cvx(data, S, D)
% Emulation of Xing's work ("Distance metric learning, with applicatin to 
% clustering with side-information) with cvx.
% Computes linear metric learning.
% 
% data: data in columns (Nxd)
% S: square NxN similarity matrix, 1 is a similar pair to be used in the
% optimization
% D: square NxN dissimilarity matrix, 1 is a dissimilar pair to be used in
% the optimization
% 
% Work protected by the Attribution-NonCommercial-ShareAlike 3.0 Unported
% Creative Commons license
% http://creativecommons.org/licenses/by-nc-sa/3.0/


% each row is a example:
s = size(data);
N = s(1);     % number of examples
d = s(2);     % dimensionality of examples            

A = eye(d,d)*0.1;
W = zeros(d,d);

%same pair groups
if size(S,1) == size(S,2) %%inefficient and very expensive for large data sets
    for i = 1:N,
        for j = i+1:N,
            if S(i,j) == 1,
                d_ij = data(i,:) - data(j,:);
                W = W + (d_ij'*d_ij);
            end
        end
    end
else %%more efficient, a list of dissimilar samples
    for sample = 1:size(S,1)
        i = S(sample,1);
        j = S(sample,2);

        d_ij = data(i,:) - data(j,:);
        W = W + (d_ij'*d_ij);
    end
end
w = unroll(W);
t = w' * unroll(A)/100;     


%different pair 
W_D = zeros(d,d);
if size(D,1) == size(D,2) %%inefficient and very expensive for large data sets
    for i = 1:N,
        for j = i+1:N,
            if D(i,j) == 1,
                d_ij = data(i,:) - data(j,:);
                W_D = W_D + (d_ij'*d_ij);
            end
        end
    end
else %%more efficient, a list of dissimilar samples
    for sample = 1:size(D,1)
        i = D(sample,1);
        j = D(sample,2);

        d_ij = data(i,:) - data(j,:);
        W_D = W_D + (d_ij'*d_ij);
    end
end
w_D = unroll(W_D);


% original Xing full matrix metric learning formulation 
cvx_begin
cvx_quiet(false);
    variable optA(d,d) symmetric;	
    optA == semidefinite(d,d); % optA >= 0; 
    maximize sum(sqrt(w_D.*unroll(optA)))
    subject to	
        sum(w'*unroll(optA)) <= t;  
cvx_end