-
Notifications
You must be signed in to change notification settings - Fork 2
/
test_greedy_algorithms.m
64 lines (64 loc) · 2.27 KB
/
test_greedy_algorithms.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
clear; close all;
%% Initial value
% m rows -> equations, n columns -> unknowns
m = 128; n = 256; sparseCardinality = 12;
a = normc(randn(m, n));
x = zeros(n, 1);
normalizedErrorBound = 1e-6;
% sparse support, or index. 'false' ensures no repeat
sparseSupport = randsample(n, sparseCardinality, 'false');
% x is the s-sparse vector
randnTemp = randn(sparseCardinality, 1);
x(sparseSupport) = randnTemp;
y = a * x;
%% Generated solutions and their normalized error
xExact = x;
xLeastSquare = a \ y;
xBasic = pinv(a) * y;
xOmp = orthogonal_matching_pursuit(sparseCardinality, a, y, normalizedErrorBound);
xSp = subspace_pursuit(sparseCardinality, a, y, normalizedErrorBound);
xIht = iterative_hardthresholding(sparseCardinality, a, y, normalizedErrorBound);
% errors
errorExact = norm(y - a * xExact) / norm(y);
errorLeastSquare = norm(y - a * xLeastSquare) / norm(y);
errorBasic = norm(y - a * xBasic) / norm(y);
errorOmp = norm(y - a * xOmp) / norm(y);
errorSp = norm(y - a * xSp) / norm(y);
errorIht = norm(y - a * xIht) / norm(y);
flag = 1;
%% Ground truth solution vs greedy solutions
xExactIndex = find (xExact ~= 0);
xOmpIndex = find(xOmp ~= 0);
xSpIndex = find(xSp ~= 0);
xIhtIndex = find(xIht ~= 0);
figure;
% orthogonal matching pursuit
ompFig = subplot(3, 1, 1);
scatter(xExactIndex, xExact(xExact ~= 0), 'x', 'black');
hold on;
xOmpScatter = scatter(xOmpIndex, xOmp(xOmp ~= 0), 'o', 'red');
xlim([0 n]);
title('Sparse solution (x): orthogonal matching pursuit');
xlabel('x basis index');
ylabel('x entries value');
legend('Exact solution', 'OMP solution', 'Location', 'bestoutside');
% subspace pursuit
spFig = subplot(3, 1, 2);
scatter(xExactIndex, xExact(xExact ~= 0), 'x', 'black');
hold on;
xSpScatter = scatter(xSpIndex, xSp(xSp ~= 0), 'o', 'magenta');
xlim([0 n]);
title('Sparse solution (x): subspace pursuit');
xlabel('x basis index');
ylabel('x entries value');
legend('Exact solution', 'SP solution', 'Location', 'bestoutside');
% iterative hardthresholding
ihtFig = subplot(3, 1, 3);
scatter(xExactIndex, xExact(xExact ~= 0), 'x', 'black');
hold on;
xIhtScatter = scatter(xIhtIndex, xIht(xIht ~= 0), 'o', 'blue');
xlim([0 n]);
title('Sparse solution (x): iterative hardthresholding');
xlabel('x basis index');
ylabel('x entries value');
legend('Exact solution', 'IHT solution', 'Location', 'bestoutside');