function lab3()

load 'MNIST_subset.mat';

%% Normalization
for i=1:size(tr_feats,1)
    tr_feats(i,:)=tr_feats(i,:)/sqrt(sum(tr_feats(i,:).^2));
end
for i=1:size(te_feats,1)
    te_feats(i,:)=te_feats(i,:)/sqrt(sum(te_feats(i,:).^2));
end

%% task 1
index2=find(tr_label==2);
index5=find(tr_label==5);
tr_label=zeros(size(tr_label));
tr_label(index2)=1;
tr_label(index5)=-1;
tr_label=tr_label([index2;index5]);
tr_feats=tr_feats([index2;index5],:);

index2=find(te_label==2);
index5=find(te_label==5);
te_label=zeros(size(te_label));
te_label(index2)=1;
te_label(index5)=-1;
te_label=te_label([index2;index5]);
te_feats=te_feats([index2;index5],:);

weights=FH_perceptron(tr_feats,tr_label);

% calculate the accuracy on the training set
R=[tr_feats ones(size(tr_label))]*weights;
Y = sign(R); Y(Y==0)=-1;
error = sum(Y~=tr_label);
fprintf('The classification accuracy on the training set is %f.\n',(size(tr_label,1)-error)/size(tr_label,1));

% calculate the accuracy on the test set
R=[te_feats ones(size(te_label))]*weights;
Y = sign(R); Y(Y==0)=-1;
error = sum(Y~=te_label);
fprintf('The classification accuracy on the test set is %f.\n',(size(te_label,1)-error)/size(te_label,1));

end










function weights=FH_perceptron(feats,label,eta)
%% feats: num_of_samples * num_of_dimension
%% eta: learning rate

if nargin<3
    eta=0.1;
end

max_iter=1000;  % maximum of iteration number, in case of non-linear data

assert(sum(abs(label)==1)==size(label,1));  % assert that label contains either 1 or -1

[nums,ndims]=size(feats);

rand('state',0);    % initialize the random seed, in order to make the rand function generate the same numbers

% initialization
weights = rand(ndims+1,1)*2-1;    % initialize weights to random numbers ranging from -1 to 1

feats = [feats ones(nums,1)];     % augment feats with another column filled with 1

delta_weights = ones(ndims+1,1);

iter=1;
while abs(sum(delta_weights))>0.1 && iter<max_iter
    R = feats*weights;
    Y = sign(R);
    Y(Y==0)=-1;
    delta_weights = eta * feats' * (label-Y);
    weights = weights + delta_weights;
%     error(iter)=abs(sum(delta_weights));
    iter=iter+1;
end

end