InverseOptimization/figure2_a_code.m at master · sorooshafiee/InverseOptimization · GitHub

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
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
clc
clearvars
rng('shuffle')

addpath('./model/')

n = 50;
m = 50;
N = 10;
N_test = 1000;
delta = 1;
radius = 1;
eta = 0;
run_count = 100;
epsilon = [(1:9)*1e-5 (1:9)*1e-4 (1:9)*1e-3 (1:10)*1e-2];
param(1:run_count) = struct('W',[],'H',[],'h',[],'C',zeros(1,m),'d',0, ...
                            'pnorm',1,'alpha',0.9,'epsilon',epsilon,   ...
                            'delta',delta,'set_theta',struct('center',[], ...
                            'radius',radius,'pnorm',inf));
data(run_count) = struct('x',[],'s',[]);
Suboptimality = zeros(length(epsilon),run_count);
Predictability = zeros(length(epsilon),run_count);

for r = 1 : run_count
    fprintf('Running iteration %d ... \n',r);

    %======================== Setting Parameters =========================%
    param(r).W = [2 * rand(m,n) - 1; eye(n); -eye(n)];
    param(r).H = [eye(m); zeros(n,m); zeros(n,m)];
    param(r).h = [zeros(m,1); -ones(n,1); -ones(n,1)];
    center = radius + 4 * radius * rand(n,1);
    sgn = randi(2,[n,1])-1;
    center(sgn == 0) = -center(sgn == 0);
    param(r).set_theta.center = center;
    theta_star = center + 2 * radius * rand(n,1) - radius;

    %========================= Generate Dataset ==========================%
    random_x = 2*rand(n,N) - 1;
    s = param(r).W(1:m,:) * random_x;
    s = s - eta;
    suboptimal = SubLinear_Model(param(r),theta_star,s);
    suboptimal_x = [suboptimal.x];

    %==================== Solve the Inverse Problem ======================%
    data(r).x   = suboptimal_x;
    data(r).s   = s;
    opt_inv     = Linear_Inverse(param(r),data(r));
    theta       = [opt_inv.theta];

    %================= Evaluate the Model on Test Data ===================%
    random_x = 2*rand(n,N_test) - 1;
    s = param(r).W(1:m,:) * random_x;
    s = s - eta;
    suboptimal = SubLinear_Model(param(r),theta_star,s);
    suboptimal_x = [suboptimal.x];
    tmp_sub     = zeros(size(theta,2),1);
    tmp_pre     = zeros(size(theta,2),1);
    for j = 1 : size(theta,2)
        opt_model   = Linear_Model(param(r),theta(:,j),s);
        tmp_sub(j)  = mean( max(theta(:,j)'* (suboptimal_x - [opt_model.x]) - delta, 0) );
        tmp_pre(j)  = mean( sqrt(sum((suboptimal_x - [opt_model.x]).^2,1)) );
    end
    Suboptimality(:,r)  = tmp_sub;
    Predictability(:,r) = tmp_pre;
end
%%
shaded = true;
fig = figure;
set(fig, 'Units', 'normalized', 'Position', [0.35, 0.25, 0.4, 0.55])
font_size = 18;

yyaxis left
h1 = semilogx(epsilon,mean(Suboptimality,2),'linewidth',4,'color',[0, 0.447, 0.741]);
set(gca, 'FontSize', font_size);
xlabel('$\varepsilon$','Interpreter','latex', 'FontSize',26);
ylabel('\delta-Suboptimality','FontSize',font_size) % left y-axis

yyaxis right
h2 = semilogx(epsilon,mean(Predictability,2),'linewidth',4,'color',[0.85, 0.325, 0.098],'linestyle','-.');
ylabel('Predictability','FontSize',font_size) % right y-axis

if shaded
    prc = 25;
    alphaa = 0.1;
    yyaxis left
    hold on
    epsilon_2 = [epsilon, flip(epsilon)];
    fill(epsilon_2,[prctile(Suboptimality,prc,2)', flip(prctile(Suboptimality,100-prc,2))'], ...
         [0, 0.447, 0.741],'LineStyle','none');
    alpha(alphaa)
    yyaxis right
    hold on
    fill(epsilon_2,[prctile(Predictability,prc,2)', flip(prctile(Predictability,100-prc,2))'], ...
         [0.85, 0.325, 0.098],'LineStyle','none');
    alpha(alphaa)
end

legend([h1 h2],{'Suboptimality','Predictability'});

cd figs
saveas(fig,'fig2-a','png')
cd ..