[Home]
L1BLDA
A Matlab code for Robust Bhattacharyya bound LDA through adaptive non-greedy algorithms. [Code]
Reference
Chun-Na Li, Yuan-Hai Shao, Zhen Wang, Nai-Yang Deng. "Robust Bhattacharyya bound linear discriminant analysis through adaptive non-greedy algorithms". Submitted, 2018.
Main Function
function [W_all] = L1BLDA(Data,FunPara)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % % % Input:
% Data.X: Data matrix. Each column vector of Data is a data point.
% Data.Y: Data label vector.
% FunPara.rho: ADMM penalty rho.
%
% % % % Eample:
% Data.X = rand(2,20);
% Data.Y = [ones(10,1);-ones(10,1)];
% FunPara.rho = 1;
% [W_all] = L1BLDA(Data,FunPara)
%
% % % % Reference:
% "Robust Bhattacharyya bound linear discriminant analysis through adaptive non-greedy algorithms".
% Chun-Na Li, Yuan-Hai Shao, Zhen Wang, Nai-Yang Deng
% Version 1.0 -- July/2018
% Written by Chun-Na Li (na1013na@163.com)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Maxstep = 20;
rho = FunPara.rho;
eps1 = 10^(-3);
eps2 = 10^(-3);
[nFea,nSmp] = size(Data.trainX);
desireDim = nFea;
W_all = cell(desireDim,1);
classLabel = unique(Data.trainY);
nClass = length(classLabel);
C_mean = zeros(nFea,nClass);
nSmpClass = zeros(1,nClass);
subX = Data.trainX;
Sw = zeros(nFea,nFea);
Sb = zeros(nFea,nFea);
Omega = 0;
% Calculate Sw
for i=1:nClass
Index(i,:) = (Data.trainY==classLabel(i));
C_mean(:,i) = mean(subX(:,Index(i,:)),2);
nSmpClass(i) = size((subX(:,Index(i,:))),2);
TempM = subX(:,Index(i,:)) - repmat(C_mean(:,i),1,nSmpClass(i));
Sw = Sw + TempM*TempM';
end
% Calculate Sb and Omega
for i = 1:nClass
for j = 1:nClass
if i eps1 || P_con2 > eps1 || P_con3 > eps1 || D_con1 > eps1 || D_con2 > eps1 || D_con3 > eps1 )
% ADMM, W
B0 = B;
for i = 1:nClass
for j = 1:nClass
if i=0) = B(B>=0) + 1/rho_new;
B(B<0) = B(B<0) - 1/rho_new;
% ADMM, Z
Z0 = Z;
for i=1:nClass
Temp = subX(:,Index(i,:));
for j=1:nSmpClass(i)
Z(i,j,:) = W'*(Temp(:,j)-C_mean(:,i))+beta(i,j);
D_con20(i,j) = norm(rho_new*(Temp(:,j)-C_mean(:,i))*(Z(i,j)-Z0(i,j))');
end
end
ka = Omega/rho_new;
Z(Z > ka) = Z(Z > ka) - ka;
Z(Z < -ka) = Z(Z < -ka) + ka;
Z(Z<=ka & Z>=-ka) = 0;
% ADMM, D
D0 = D;
D = W - Gamma;
Po_con1 = P_con1; Po_con2 = P_con2; Po_con3 = P_con3;
Do_con1 = D_con1; Do_con2 = D_con2; Do_con3 = D_con3;
for i=1:nClass
xi = subX(:,Index(i,:));
for j=1:nSmpClass(i)
tempijbeta = W'*(xi(:,j)-C_mean(:,i))-Z(i,j);
beta(i,j,:) = beta(i,j)+tempijbeta;
P_con20(i,j) = norm(tempijbeta);
end
end
Gamma = Gamma + (D-W);
P_con1 = max(max(P_con10));
P_con2 = max(max(P_con20));
P_con3 = norm(D-W);
D_con1 = max(max(D_con10));
D_con2 = max(max(D_con20));
D_con3 = norm(rho_new*(D-D0));
P_con = max(P_con1,max(P_con2,P_con3));
D_con = max(D_con1,max(D_con2,D_con3));
mu=10;t_incr=2;t_decr=2;
if P_con>mu*D_con
rho_new=t_incr*rho_new;
elseif D_con>mu*P_con
rho_new=rho_new/t_decr;
else
rho_new=rho_new;
end
if abs(Po_con1 - P_con1) < eps2 && abs(Po_con2 - P_con2) < eps2 && abs(Po_con3 - P_con3) < eps2 < eps2 && abs(Do_con1 - D_con1) < eps2 && abs(Do_con2 - D_con2) < eps2 && abs(Do_con3 - D_con3) < eps2
break;
end
if step>Maxstep
break
end
step = step +1;
TP_con1(step) = P_con1;
TP_con2(step) = P_con2;
TP_con3(step) = P_con3;
TD_con1(step) = D_con1;
TD_con2(step) = D_con2;
TD_con3(step) = D_con3;
end
W_all{d} = W;
clear W B alpha Z beta D Gamma
else % when d = n
for i = 1:nClass
for j = 1:nClass
if i eps1 || P_con2 > eps1 || P_con3 > eps1 || D_con1 > eps1 || D_con2 > eps1 || D_con3 > eps1 )
% ADMM, W
B0 = B;
for i = 1:nClass
for j = 1:nClass
if i=0) = B(B>=0) + 1/rho_new;
B(B<0) = B(B<0) - 1/rho_new;
% ADMM, Z
Z0 = Z;
for i=1:nClass
Temp = subX(:,Index(i,:));
for j=1:nSmpClass(i)
Z(i,j,:) = W'*(Temp(:,j)-C_mean(:,i))+beta(i,j);
D_con20(i,j) = norm(rho_new*(Temp(:,j)-C_mean(:,i))*(Z(i,j)-Z0(i,j)));
end
end
ka = Omega/rho_new;
Z(Z > ka) = Z(Z > ka) - ka;
Z(Z < -ka) = Z(Z < -ka) + ka;
Z(Z<=ka & Z>=-ka) = 0;
% ADMM, D
D0 = D;
D = W - Gamma;
Po_con1 = P_con1; Po_con2 = P_con2; Po_con3 = P_con3;
Do_con1 = D_con1; Do_con2 = D_con2; Do_con3 = D_con3;
for i=1:nClass
xi = subX(:,Index(i,:));
for j=1:nSmpClass(i)
tempijbeta = W'*(xi(:,j)-C_mean(:,i))-Z(i,j);
beta(i,j,:) = beta(i,j)+tempijbeta;
P_con20(i,j) = norm(tempijbeta);
end
end
Gamma = Gamma + (D-W);
P_con1 = max(max(P_con10));
P_con2 = max(max(P_con20));
P_con3 = norm(D-W);
D_con1 = max(max(D_con10));
D_con2 = max(max(D_con20));
D_con3 = norm(rho_new*(D-D0));
P_con = max(P_con1,max(P_con2,P_con3));
D_con = max(D_con1,max(D_con2,D_con3));
mu=10;t_incr=2;t_decr=2;
if P_con>mu*D_con
rho_new=t_incr*rho_new;
elseif D_con>mu*P_con
rho_new=rho_new/t_decr;
else
rho_new=rho_new;
end
if abs(Po_con1 - P_con1) < eps2 && abs(Po_con2 - P_con2) < eps2 && abs(Po_con3 - P_con3) < eps2 < eps2 && abs(Do_con1 - D_con1) < eps2 && abs(Do_con2 - D_con2) < eps2 && abs(Do_con3 - D_con3) < eps2
break;
end
if step>Maxstep
break
end
step = step +1;
TP_con1(step) = P_con1;
TP_con2(step) = P_con2;
TP_con3(step) = P_con3;
TD_con1(step) = D_con1;
TD_con2(step) = D_con2;
TD_con3(step) = D_con3;
end
W_all{d} = W;
end
end
function [W] = QPSM(G,A)
itmax = 10;
eps = 10^-4;
it2 = 0;
[m,k] = size(A);
[~,D_G] = eig(G);
D_G = diag(D_G);
alpha = max(D_G) + 0.01;
W = [eye(k);zeros(m-k,k)];
while it2 < itmax
it2 = it2 + 1;
W0 = W;
M = 2*(alpha*eye(m) - G)*W - 2*A;
[U,~,V] = svd(M,0);
W = U*V';
if norm(W - W0) < eps
break;
end
end
Any question or advice please email to na1013na@163.com or shaoyuanhai21@163.com.
- Last updated: July 24, 2018