Commit dcaf0995 authored by Elvira Clement's avatar Elvira Clement
Browse files

Updating experiment 1

parent bcbebe3d
python xp_a_accuracy_sol.py --id SIAM
python xp_b_screening.py --id SIAM
python xp_c_viz_paper.py --id SIAM
\ No newline at end of file
# -*- coding: utf-8 -*-
import numpy as np
from xps.SIAM.setup import Setup
from xps.SIAM.xp_1_balls import xpparams
def process(setup, log=True):
# ---- load ----
solutions_filename = f"results/xp_setup{setup.setup_id}.npz"
solutions = np.load(solutions_filename, allow_pickle=True)
screenings_filename = f"results/xp_setup{setup.setup_id}_screening.npz"
screenings = np.load(screenings_filename, allow_pickle=True)
# mat_results_gap[1, i_dic, i_seq, rep, i_ratio]
mat_pvopt = solutions["mat_pvopt"]
mat_nb_zero_detected = screenings['mat_nb_zero_detected']
nb_test = mat_nb_zero_detected.shape[0]
list_tests = screenings["list_test"]
# ---- log ----
if log:
print("Experiment info")
# print(f"- run with version {solutions["version"]}")
print(f"- setup{setup.setup_id}")
print("")
# ---- processing ----
mat_pc_detected =np.zeros(
(nb_test, xpparams.nb_point, setup.nb_dic, setup.nb_sequence, setup.nb_ratio_lbd, setup.n_rep)
)
for i_dic in range(setup.nb_dic):
for i_seq in range(setup.nb_sequence):
for i_ratio in range(setup.nb_ratio_lbd):
for rep in range(setup.n_rep):
nb_0 = np.sum(mat_pvopt[i_dic, i_seq, i_ratio, rep, :] == 0)
if nb_0 > 0:
mat_detected = mat_nb_zero_detected[:, :, i_dic, i_seq, i_ratio, rep]
mat_pc_detected[:, :, i_dic, i_seq, i_ratio, rep] = mat_detected / float(nb_0)
# ---- return ----
return {
"mat_pc_detected": np.mean(mat_pc_detected, axis=5),
"list_tests": list_tests,
}
# {
# "results_mat_nb": results_mat_nb,
# "results_mat_nbz": results_mat_nbz,
# }
\ No newline at end of file
# -*- coding: utf-8 -*-
import time, argparse, sys
import numpy as np
# Algorithm import
from src import __version__
from src.solver.parameters import SlopeParameters, EnumLipchitzOptions
from src.solver.slope import primal_func, dual_func, slope_gp
# Generative models
from src.utils import get_lambda_max, gamma_sequence_generator
from src.dictionaries import generate_dic
# Screening
from src.screening.gap_test_all import GapTestAll
# XP import
from experiments.SIAM.slopepb import SlopePb
from experiments.SIAM.setup import Setup
from experiments.SIAM.xp_1_balls import xpparams
parser=argparse.ArgumentParser()
parser.add_argument('--erase', help='erase existing results', action="store_true")
parser.add_argument('--id', help='setup id', type=str)
args=parser.parse_args()
# -------------------------
# Load Setup
# -------------------------
setup = Setup(args.id)
out_file_name = f"results/xp_setup{args.id}.npz"
mat_seed = np.random.randint(
0, 2**32-1,
size=(setup.nb_dic, setup.nb_sequence, setup.n_rep),
)
mat_pvopt = np.zeros(
(setup.nb_dic, setup.nb_sequence, setup.nb_ratio_lbd, setup.n_rep, setup.n)
)
mat_dvopt = np.zeros(
(setup.nb_dic, setup.nb_sequence, setup.nb_ratio_lbd, setup.n_rep, setup.m)
)
try:
if args.erase:
raise FileNotFoundError
# Try to load existing results
load_results = np.load(out_file_name, allow_pickle=True)
mat_seed = load_results["mat_seed"]
mat_pvopt = load_results["mat_pvopt"]
mat_dvopt = load_results["mat_dvopt"]
except FileNotFoundError:
pass
# -------------------------
# Xp
# -------------------------
# For each trial
# 1. Compute high accuracy solution
# 2. Create ideal safe ball
# 3. increase radius
nb_xp = setup.nb_dic * setup.n_rep * setup.nb_sequence * setup.nb_ratio_lbd
i_xp = 0
for i_dic in range(setup.nb_dic):
for i_seq, seq in enumerate(setup.list_sequence):
for rep in range(setup.n_rep):
np.random.seed(mat_seed[i_dic, i_seq, rep])
# ---- 1. Gen data ----
matA = generate_dic(
setup.list_dic[i_dic],
setup.m,
setup.n,
setup.normalize
)
lip = np.linalg.norm(matA, ord=2)**2
vecy = np.random.randn(setup.m)
vecy /= np.linalg.norm(vecy, 2)
# ---- 2. Compute parameters ----
vec_gammas = gamma_sequence_generator(
setup.m,
setup.n,
setup.list_sequence[i_seq][0],
setup.list_sequence[i_seq][1:]
)
lbd_max = get_lambda_max(vecy, matA, vec_gammas)
# ---- 3. XP ----
for i_ratio, ratio in enumerate(setup.list_ratio_lbd):
i_xp += 1
slopePb = SlopePb(matA, vecy, vec_gammas, ratio, lbdmax=lbd_max)
vecx_hat = mat_pvopt[i_dic, i_seq, i_ratio, rep, :]
vecu_hat = mat_dvopt[i_dic, i_seq, i_ratio, rep, :]
gap = slopePb.eval_gap(vecx_hat, vecu_hat)
print(f"xp {i_xp} / {nb_xp} --- (gap={gap})")
if (gap > xpparams.stopping_gap):
gap_old = gap
# ---- 3a. Find solution ----
params = SlopeParameters()
params.vecx_init = np.copy(vecx_hat)
params.lipchitz_constant = lip
params.lipchitz_update = EnumLipchitzOptions.EXACT
params.max_it = 1e7
params.gap_stopping = xpparams.stopping_gap
params.time_stopping = np.inf
params.screening1 = GapTestAll(vec_gammas)
params.eval_gap = True
params.eval_gap_it = setup.eval_gap_it
params.accelerated = False
params.verbose = False
out_slope = slope_gp(vecy, matA, ratio * lbd_max, vec_gammas, params)
dv = vecy - matA @ out_slope["sol"]
dv = slopePb.make_dual_scaling(dv)
gap = slopePb.eval_gap(out_slope["sol"], dv)
if gap <= gap_old:
mat_pvopt[i_dic, i_seq, i_ratio, rep, :] = out_slope["sol"]
mat_dvopt[i_dic, i_seq, i_ratio, rep, :] = dv
# Save
np.savez(out_file_name,
mat_seed=mat_seed,
mat_pvopt=mat_pvopt,
mat_dvopt=mat_dvopt,
version = __version__,
allow_pickle=True
)
# if i_xp ==3:
# exit()
\ No newline at end of file
# -*- coding: utf-8 -*-
import argparse
import numpy as np
# Algorithm import
from src import __version__
from src.dictionaries import generate_dic
from src.utils import get_lambda_max, gamma_sequence_generator
# Screening
from src.screening.gap_test_p_1 import GapTestPequalOne
from src.screening.gap_test_p_q import GapTestPequalQ
from src.screening.gap_test_all import GapTestAll
# XP import
from experiments.SIAM.slopepb import SlopePb
from experiments.SIAM.setup import Setup
from experiments.SIAM.xp_1_balls import xpparams
parser=argparse.ArgumentParser()
parser.add_argument('--id', help='setup id', type=str, default=1)
args=parser.parse_args()
# -------------------------
# Load Setup
# -------------------------
setup = Setup(args.id)
solutions_filename = f"results/xp_setup{args.id}.npz"
try:
# Try to load existing results
load_results = np.load(solutions_filename, allow_pickle=True)
mat_seed = load_results["mat_seed"]
mat_pvopt = load_results["mat_pvopt"]
mat_dvopt = load_results["mat_dvopt"]
except FileNotFoundError:
print("No result found")
sys.exit(1)
# -------------------------
# Screening Xp
# -------------------------
# percentage of detected zero
mat_nb_zero_detected = np.full(
(3, xpparams.nb_point, setup.nb_dic, setup.nb_sequence, setup.nb_ratio_lbd, setup.n_rep),
np.nan,
dtype=int
)
screening_filename = f"results/xp_setup{args.id}_screening.npz"
i_xp = 0
nb_xp = setup.nb_dic * setup.n_rep * setup.nb_sequence * setup.nb_ratio_lbd
for i_dic in range(setup.nb_dic):
for i_seq, seq in enumerate(setup.list_sequence):
for rep in range(setup.n_rep):
np.random.seed(mat_seed[i_dic, i_seq, rep])
# ---- 1. Gen data ----
matA = generate_dic(
setup.list_dic[i_dic],
setup.m,
setup.n,
setup.normalize
)
vecy = np.random.randn(setup.m)
vecy /= np.linalg.norm(vecy)
# ---- 2. Compute parameters ----
vec_gammas = gamma_sequence_generator(
setup.m,
setup.n,
setup.list_sequence[i_seq][0],
setup.list_sequence[i_seq][1:]
)
lbd_max = get_lambda_max(vecy, matA, vec_gammas)
# ---- 3. XP ----
for i_ratio, ratio in enumerate(setup.list_ratio_lbd):
i_xp += 1
print(f"xp {i_xp} / {nb_xp}")
slopePb = SlopePb(matA, vecy, vec_gammas, ratio, lbdmax=lbd_max)
vecx_hat = mat_pvopt[i_dic, i_seq, i_ratio, rep, :]
vecu_hat = mat_dvopt[i_dic, i_seq, i_ratio, rep, :]
gap = slopePb.eval_gap(vecx_hat, vecu_hat)
Atu = matA.T @ vecu_hat
# ---- 3b. Build thin safe ball ----
rgap = np.sqrt(2 * gap)
# ---- 3c. Testing sphere ----
list_tests = [
GapTestPequalQ(),
GapTestAll(vec_gammas),
GapTestPequalOne(vec_gammas, np.arange(setup.n, dtype=np.double)),
]
for i_offset, offset in enumerate(xpparams.vec_offsets):
for i_test, test in enumerate(list_tests):
out = test.apply_test(np.abs(Atu), gap, ratio * lbd_max, vec_gammas, offset_radius=offset)
mat_nb_zero_detected[i_test, i_offset, i_dic, i_seq, i_ratio, rep] = np.sum(out)
# Save
np.savez(
screening_filename,
mat_nb_zero_detected=mat_nb_zero_detected,
list_test = [test.get_name() for test in list_tests],
version = __version__,
allow_pickle=True
)
\ No newline at end of file
# -*- coding: utf-8 -*-
from decimal import Decimal
import json, argparse
import numpy as np
# XP import
from xps.SIAM.slopepb import SlopePb
from xps.SIAM.setup import Setup
from xps.SIAM.xp_1_balls import xpparams
from xps.SIAM.xp_1_balls.process_data import process
parser=argparse.ArgumentParser()
parser.add_argument('--noshow', help='do not display figures', action="store_true")
parser.add_argument('--save', help='save figure', action="store_true")
parser.add_argument('--id', help='setup id', type=str, default="SIAM")
parser.add_argument('--i_seq', type=int, default=0)
args=parser.parse_args()
import matplotlib
if args.noshow:
matplotlib.use('PS')
else:
matplotlib.use("TkAgg")
import matplotlib.pyplot as plt
from matplotlib.legend_handler import HandlerBase
import matplotlib.font_manager as font_manager
# -------------------------
# Font stuff
# -------------------------
fs = 20
font_math = font_manager.FontProperties(
fname='../fonts/cmunrm.ttf',
# weight='normal',
# style='normal',
math_fontfamily="cm",
size=fs+2
)
font_text = font_manager.FontProperties(
fname='../fonts/cmunrm.ttf',
# weight='normal',
# style='normal',
math_fontfamily="cm",
size=fs+2
)
font_ttt = font_manager.FontProperties(
# fname='../fonts/ectt1000.ttf',
fname='../fonts/computer-modern/cmuntt.ttf',
weight='bold',
style='normal',
size=fs
)
# -------------------------
# Load Results
# -------------------------
setup_oscar = Setup(args.id)
dic_process_oscar = process(setup_oscar)
mat_pc_detected_oscar = dic_process_oscar["mat_pc_detected"]
list_tests_oscar = dic_process_oscar["list_tests"]
# -------------------------
# Plot Results
# -------------------------
i_dic = 2
i_lbd = 1
fs=22
fs_ylabels = 20
list_colors = ["tab:blue", "tab:orange", "tab:green"]
list_legends = ["$p_q=q\\;\\forall q$", "all", "$p_q=1\\;\\forall q$"]
# "best $r_q \\;\\forall q$ "
print("printing xp_0_ball parameters with")
print(f"- {setup_oscar.list_dic[i_dic]} dictionary")
print(f"- lbd / lbd_max = {setup_oscar.list_ratio_lbd[i_lbd]}")
f, ax = plt.subplots(1, 1, figsize=(.7*16, .6*9), sharex=True, sharey=True)
# ax.set_title(f"OSCAR-{i_seq+1}", fontsize=fs+2)
ax.set_xlabel(
r"$R$",
fontsize=fs+2,
fontproperties=font_math,
)
ax.set_ylabel(
"% of zero entries detected",
fontsize=fs+2,
fontproperties=font_text
)
for tick in ax.xaxis.get_major_ticks():
tick.label.set_fontproperties(font_math)
tick.label.set_fontsize(20)
for tick in ax.yaxis.get_major_ticks():
tick.label.set_fontproperties(font_math)
tick.label.set_fontsize(20)
for i_test in [2, 0, 1]:
ax.plot(
xpparams.vec_offsets,
100 * mat_pc_detected_oscar[i_test, :, i_dic, args.i_seq, i_lbd],
label = list_legends[i_test],
linewidth=4.,
alpha=.9,
color=list_colors[i_test]
)
ax.legend(
fontsize=fs-2,
prop=font_math
)
ax.set_xscale("log")
ax.set_xlim([1e-6, 1e0])
ax.set_ylim([-2, 102])
if args.save:
filename = f"figs/xp_illustration_screening{args.i_seq}.eps"
plt.savefig(filename, bbox_inches='tight')
if not args.noshow:
plt.show()
\ No newline at end of file
python xp_a_accuracy_sol.py --id 1a
python xp_b_screening.py --id 1a
python xp_c_viz_fig1.py --noshow --save
python xp_c_viz_fig2.py --noshow --save
\ No newline at end of file
python xp_a_accuracy_sol.py --id SIAM
python xp_b_screening.py --id SIAM
python xp_c_viz_fig1.py --id SIAM --save
python xp_c_viz_fig2.py --id SIAM --save
\ No newline at end of file
# -*- coding: utf-8 -*-
import time, argparse, sys
from pathlib import Path
import numpy as np
......@@ -33,6 +34,7 @@ args=parser.parse_args()
# -------------------------
setup = Setup(args.id)
Path("results").mkdir(parents=True, exist_ok=True)
out_file_name = f"results/xp_setup{args.id}.npz"
mat_seed = np.random.randint(
......
# -*- coding: utf-8 -*-
import argparse
from pathlib import Path
import numpy as np
......@@ -132,7 +133,10 @@ ax.set_xlim([1e-6, 1e0])
ax.set_ylim([-2, 102])
if args.save:
filename = f"figs/xp_illustration_screening{args.i_seq}.eps"
folderfig = "figs"
Path(folderfig).mkdir(parents=True, exist_ok=True)
filename = folderfig + f"/xp_illustration_screening{args.i_seq}.eps"
plt.savefig(filename, bbox_inches='tight')
if not args.noshow:
......
# -*- coding: utf-8 -*-
import argparse
from pathlib import Path
import numpy as np
......@@ -170,6 +171,9 @@ for i_seq in range(3):
# exit()
if args.save:
filename = f"figs/xp0_{setup_oscar.list_dic[i_dic]}.eps"
folderfig = "figs"
Path(folderfig).mkdir(parents=True, exist_ok=True)
filename = folderfig + f"/xp0_{setup_oscar.list_dic[i_dic]}.eps"
# plt.rcParams['pdf.fonttype'] = 42
plt.savefig(filename, bbox_inches='tight')
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment