<jupyter_start><jupyter_text>Figure 3(i)This notebook contains three major parts: - I. initialize test instances based on baseline lps generated from previous notebook- II. Experiment (lean c using deep_inv_opt package)- III. Data Analysis (generating boxplot as Figure 3(i) in our paper)<jupyter_code># Copyright (C) to <NAME>, <NAME>, <NAME>. All Rights Reserved.

import numpy as np
import scipy.optimize as opt
from scipy.spatial import ConvexHull
import matplotlib.pyplot as plt
import deep_inv_opt as io
import os
from shutil import copyfile
import time
import torch

def seed_random(seed=0):
    np.random.seed(seed)

# All generated files will be saved in the following directory
PATH = "C:/Users/Public/Downloads/"<jupyter_output><empty_output><jupyter_text>I. initialize vertex taget for task of learning cBased on instance generated from **LP generator**, we sample initial c vectors (c_init) and vertex targets.- c_init is sampled through normal distribution, c_i = Norm(0,1)- vertex targets are randomly selected from all optimal solutions from the LP instances generated above.<jupyter_code>def read_baselineLP(filename):
    lines=[]
    with open(filename, "rt") as fp:
        for line in fp:
            lines.append(line.rstrip('\n')) 
    
    m,n = np.fromstring(lines[0], dtype=int, sep=' ')
    
    temp = [lines[i+1: i+1+m] for i in range(len(lines)) if "A_ub" in lines[i]]
    A_ub=np.zeros((m,n))
    b_ub=np.zeros((m,1))
    for j in range(m):
        A_ub[j] = np.fromstring(temp[0][j], dtype=float, sep=' ')
    
    temp = [lines[i+1: i+1+m] for i in range(len(lines)) if "b_ub" in lines[i]]
    for j in range(m):
        b_ub[j] = np.fromstring(temp[0][j], dtype=float, sep=' ')
        
    temp = [lines[i] for i in range(len(lines)) if "vertices" in lines[i]]
    num_vertices = int(temp[0].split()[0])

    temp = [lines[i+1: i+1+num_vertices] for i in range(len(lines)) if "vertices" in lines[i]]
    vertices = np.zeros((num_vertices,n))
    for j in range(num_vertices):
        vertices[j] = np.fromstring(temp[0][j], dtype=float, sep=' ')

    return A_ub, b_ub, vertices


""" uncomment the following lines to run a quick test of this cell  """

# print("===== Testing \"read_BaselineLP\"=====")
# nVar, nCons, ind_lp = 2,4,0
# LPbaseline = PATH + "Deep_Inv_Opt/LP baseline/n=%d,m=%d/%svar_LP%s.txt"%(nVar, nCons, str(nVar).zfill(2),str(ind_lp+1).zfill(2))
# A_ub, b_ub, vertices = read_baselineLP(LPbaseline)
# print(A_ub, "\n", b_ub, "\n", vertices)
def Initialize_3a_instance(ind_lp, nCons, nVar):
    def sample_c (num, n, dist_type):
        if dist_type == "uniform":
            temp = np.random.uniform(-1,1,(num,n))
        elif dist_type == "normal":
            mu = 0.0
            sigma = 1.0
            temp = np.random.normal(mu, sigma, (num,n))
        return temp
    
    LPbaseline = PATH + "Deep_Inv_Opt/LP baseline/n=%d,m=%d/%svar_LP%s.txt"%(nVar, nCons, str(nVar).zfill(2),str(ind_lp+1).zfill(2))
    _, _, opt_sol = read_baselineLP(LPbaseline)
    vertex_ind = np.random.randint(0,len(opt_sol)) #randomly select target
    c_init = sample_c(1, nVar, "normal")

    directory = PATH + "Deep_Inv_Opt/figure3_i/n=%d,m=%d/"%(nVar, nCons)
    if not os.path.exists(directory):
        os.makedirs(directory)
    LP3a = directory + "%svar_LP%s.txt"%(str(nVar).zfill(2), str(ind_lp+1).zfill(2))

    copyfile(LPbaseline, LP3a)
    with open(LP3a, "a") as file:
        # print LP size, m - # constrains, n - # # var
        file.write("vertex target %d\n"%(vertex_ind+1))
        np.savetxt(file, opt_sol[vertex_ind].reshape(1, nVar), fmt='%.6f')
        file.write("\nc_init\n")
        np.savetxt(file, c_init.reshape(1, nVar), fmt='%.6f')
        file.write("\n")    

""" uncomment the following lines to run a quick test of this cell  """
# print("===== Testing \"Initialize_3a_instance\"=====")
# nVar, nCons, ind_lp = 2,4,0
# Initialize_3a_instance(ind_lp, nCons, nVar)<jupyter_output><empty_output><jupyter_text>Run the following cell to intialize instances for Figure3(i)<jupyter_code>seed_random()

_RANGE_CONS = [[4,  8,  16],
               [20, 36, 80]]
_RANGE_VAR = [2, 10]

NUM_INS = 50

for i, nVar in enumerate(_RANGE_VAR):
    for nCons in _RANGE_CONS[i]:
        for ind_lp in range(NUM_INS):
            Initialize_3a_instance(ind_lp, nCons, nVar)<jupyter_output><empty_output><jupyter_text>II. Experiment -  Learn c using Deep_Inv_Opt package Utility Functions- **ADG**: calcualte Error (Abs duality gap)- **MSE**: calculate Error (Square Error)- **Read_3aLP**: read test instance for experiment/data analysis<jupyter_code>def ADG(c,x_predict, x_target):
    c_vector = io.tensor(c).view(-1,1)
    x_predict = io.tensor(x_predict).view(-1,1)
    x_target = io.tensor(x_target).view(-1,1)
    err = torch.abs(c_vector.t() @ (x_target - x_predict))
    
    return err.detach().numpy().ravel()

def MSE(x1, x2):
    x1 = io.tensor([x1]).view(-1,1)
    x2 = io.tensor([x2]).view(-1,1)
    err = torch.mean(torch.sum((x1 - x2)**2, 0))
    
    return err.detach().numpy().ravel()

def Read_3aLP(filename):    
    A_ub, b_ub, _ = read_baselineLP(filename)
        
    lines=[]    
    with open(filename, "rt") as fp:
        lines=fp.readlines()
        for line in fp:
            lines.append(line.rstrip('\n')) 

    temp = [lines[j+1] for j in range(len(lines)) if "c_init" in lines[j]]
    c_init = np.fromstring(temp[0], dtype=float, sep=' ') 

    temp = [lines[j+1] for j in range(len(lines)) if "vertex target" in lines[j]]
    x_target = np.fromstring(temp[0], dtype=float, sep=' ') 

    return (A_ub, b_ub, x_target, c_init)

""" uncomment the following lines to run a quick test of this cell  """
# print("===== Testing \"Read_3aLP\"=====")
# nVar, nCons, ind_lp = 10,20,0
# filename = PATH + "Deep_Inv_Opt/figure3_i/n=%d,m=%d/%svar_LP%s.txt"%(nVar, nCons, str(nVar).zfill(2),str(ind_lp+1).zfill(2))
# print(Read_3aLP(filename))<jupyter_output><empty_output><jupyter_text>Random Hyperparameter SearchSample values of hyperparmeters from pre-defined sets, which are to be used in random hyper-parameter search<jupyter_code># pre-defined set of hyperparameters
_RANGE_T0 = [0.5, 1, 5, 10]
_RANGE_MU= [1.5, 2, 5, 10, 20]
_RANGE_LR_C = [1,10 ,100, 1000]
_RANGE_EPS = [1e-5]

def sample_hyperparam(filename, num_set = 20, print_hyperparam = True):
    hyperParam = [[]]*num_set
    for i in range(20):
        mu = _RANGE_MU[np.random.randint(0, len(_RANGE_MU))]
        t0 = _RANGE_T0[np.random.randint(0, len(_RANGE_T0))]
        eps = _RANGE_EPS[np.random.randint(0, len(_RANGE_EPS))]
        lr_c = _RANGE_LR_C[np.random.randint(0, len(_RANGE_LR_C))]
        hyperParam[i] = (mu, t0, eps, lr_c)
        if print_hyperparam:
            print("hyperparam set %d"%(i))
            print(hyperParam[i])
        with open(filename, 'w') as file:
            for j in range(len(hyperParam)):
                file.write("hyperParam%s: "%(str(j+1).zfill(2)))
                np.savetxt(file, np.mat(hyperParam[j]).ravel(), fmt='%.6f')
            file.write("\n")

def read_hyperParam(filename, num_set = 20):
    with open(filename, 'rt') as file:
        lines=file.readlines()
        temp = [lines[i] for i in range(len(lines)) if "hyperParam" in lines[i]]
        hyperParam = [[]]*num_set
        for i in range(len(temp)):
            te = temp[i]
            te = te[14:]
            te = te[:-2]
            te = te.split(' ')
            hyperParam[i] = [float(j) for j in te]

    return hyperParam

print("==== Sample Hyperparameters values for Random Hyperparameter Search====")
seed_random()
filename = PATH + "Deep_Inv_Opt/figure3_i/Hyper_Param.txt"
print_hyperparam = True # set this to False to disable the printing
sample_hyperparam(filename)
read_hyperParam(filename)<jupyter_output><empty_output><jupyter_text>Learn with Deep_Inv_Opt package using random hyperparameter search Please Note, finishing the entire experiments (i.e., 6 * 50 instances for Training withboth of ADG and MSE) might take 1-2 days. You can see our paper for the experiment results directly.<jupyter_code>def exp_randomsearch(ind_lp, LPInstance, hyperParam):
    #function to record result of each hyperparameter search
    def record_result(filename, hyperpara, result):
        loss, time, c_final, x_final = result
        with open(filename, 'a') as file:
            mu, t0, eps, lr_c = hyperParam[ind_search]
            temp = (mu, t0, eps, lr_c, loss, time)
            string = ("hyperparam"+str(ind_search+1).zfill(2)+" "+ str(temp) )
            file.write(str(string))
            file.write("\n")
            file.write("c_final ")
            np.savetxt(file, c.reshape(1,nVar), fmt='%.6f') 
            file.write("x_final ")
            np.savetxt(file, x.reshape(1,nVar), fmt='%.6f')     

    A_ub, b_ub, x_target, c_init = LPInstance
    nCons,nVar = A_ub.shape
    A_ub = io.tensor(A_ub)
    b_ub = io.tensor(b_ub)  
    c_init = io.tensor(c_init).view(-1,1)
    x_target = io.tensor(x_target).view(-1,1)

    print("||||| Ins %d (%d var, %d constraint), vertex target: "
          %(ind_lp+1, nVar, nCons))
    print("||||| vertex target ",x_target.detach().numpy())
    
    for loss in ["ADG", "MSE"]:   
        # record LP characteristics
        _RESULT_FILE = PATH + 'Deep_Inv_Opt/figure3_i/EXP_%svar_%scon_%s.txt'%(str(nVar).zfill(2),str(nCons).zfill(2), loss)
        with open(_RESULT_FILE, 'a') as file:
            file.write( "LP%s, %svar %scon:\n"
                       %(str(ind_lp+1).zfill(2),str(nVar).zfill(2),str(nCons).zfill(2)))
            file.write( "Vertex \n" )
            np.savetxt(file, x_target.detach().numpy().reshape(1,nVar), fmt='%.6f')
            file.write("c_init: ")
            np.savetxt(file, c_init.numpy().reshape(1,nVar), fmt='%.6f')                       

        # initialize loss function
        if loss =="ADG":
            loss_callback = io.abs_duality_gap
            eps_decay = False

        elif loss =="MSE":
            loss_callback = io.squared_error 
            eps_decay = True
        
        for ind_search in range(len(hyperParam)):
            # initialize hyper parameters
            mu, t0, eps, lr_c = hyperParam[ind_search]
            print("loss func %s, HyperParam %d = (%.1f, %.1f, %.5f, %.1f"
                  %(loss, ind_search+1, mu, t0, eps, lr_c),")")
            print("---- c_init ", c_init.numpy().ravel())
            print("---- Hyper-Param Search %d / %d: "%(ind_search+1, len(hyperParam)))
            
            solver = io.custom_linprog(t0 = t0, mu = mu, max_steps = 10000)

            start_time = time.time()
            try:
                c_final, _, _, err, x_predict, init_err = io.inverse_linprog(x_target, 
                                                                            c_init,  A_ub, b_ub,
                                                                            max_steps=200,
                                                                            eps = eps, #defaul value for the end of eps_decay
                                                                            eps_decay = eps_decay,
                                                                            learn_rate_c=lr_c,
                                                                            loss=loss_callback,
                                                                            return_loss=True,
                                                                            solver=solver)
                runtime = time.time() - start_time

                c = c_final.detach().numpy().ravel()
                x = x_predict.detach().numpy().ravel()
                l = err.detach().numpy().ravel()
                t = runtime
                
                print("---- Final Loss", l)

                #record Hyper param values
                result = (l[0], t, c, x)
                record_result(_RESULT_FILE, hyperParam[ind_search], result)   
                
            except: # catch error and return failed hyperparam search
                print("find an error, return failed hyperparam search")
                c = np.ones((c_init.shape))*999
                x = np.ones((c_init.shape))*999              
            
                result = (999, 999, c, x)
                record_result(_RESULT_FILE, hyperParam[ind_search], result) 
                pass
            print()

            
RANGE_CONS = [[4,  8,  16],
               [20, 36, 80]]
RANGE_VAR = [2, 10]

NUM_INS = 50

for i, nVar in enumerate(RANGE_VAR):
    for nCons in RANGE_CONS[i]:
        for ind_lp in range(NUM_INS):
            filename = PATH + "Deep_Inv_Opt/figure3_i/n=%d,m=%d/%svar_LP%s.txt"%(nVar, nCons, str(nVar).zfill(2),str(ind_lp+1).zfill(2))
            LPInstance = Read_3aLP(filename)
            hyperParam = read_hyperParam(PATH + 'Deep_Inv_Opt/figure3_i/Hyper_Param.txt')
            exp_randomsearch(ind_lp, LPInstance, hyperParam)<jupyter_output><empty_output><jupyter_text>III. Data Analysis Unility Function- **read_finalresult**: read result from data file<jupyter_code>def read_finalresult(num_lp, nCons, nVar, loss):
    print("\tloss = ", loss)
    c_init = np.zeros((num_lp, nVar))

    # read c_init and LP ins
    for ind_lp in range(num_lp):
        filename = PATH + "Deep_Inv_Opt/figure3_i/n=%d,m=%d/%svar_LP%s.txt"%(nVar, nCons, str(nVar).zfill(2),str(ind_lp+1).zfill(2))
        A, b, opt_sol, c = Read_3aLP(filename)
        c_init[ind_lp] = c

    # read x_target
    x_target = np.zeros((num_lp, nVar))
    result = []
    c_final = []
    x_final = []

    with open(PATH + 'Deep_Inv_Opt/figure3_i/FinalExp _%svar_%scon_%s.txt'
              %(str(nVar).zfill(2),str(nCons).zfill(2), loss), 'rt') as file:
        lines=file.readlines()
        temp = [lines[i+1] for i in range(len(lines)) if "Vertex" in lines[i]]
        for i in range((len(temp))):
            te =temp[i]
            te = te.split(' ')
            te = [float(d) for d in te]
            x_target[i] = te
        #read exp results of every trail (mu, to, eps, lr_c, err, runtime)
        temp = [lines[i] for i in range(len(lines)) if "hyperparam" in lines[i]]
        for i in range(len(temp)):
            te = temp[i]
            te = te[14:]
            te = te[:-2]

            te = te.split(',')
            te = [float(d) for d in te]
            result.append(te)
        result= np.mat(result)
        #from result find the index of the best trail for each LP ins
        best_ind = []
        for i in range(num_lp):
            best_ind.append(int ( i*20 + np.argmin(result[i*20:(i+1)*20,-2])))
        best_ind = np.array(best_ind)
        # read c_final for every trail
        temp = [lines[i] for i in range(len(lines)) if "c_final" in lines[i]]
        for i in range(len(temp)):
            te = temp[i]
            te = te[8:]                
            te = te.split(' ')
            te = [float(d) for d in te]
            c_final.append(te)
        c_final= np.mat(c_final)
        
        # read x_final for every trail
        temp = [lines[i] for i in range(len(lines)) if "x_final" in lines[i]]
        for i in range(len(temp)):
            te = temp[i]
            te = te[8:]                
            te = te.split(' ')
            te = [float(d) for d in te]
            x_final.append(te)
        x_final= np.mat(x_final)
        
    return result, c_init, c_final, best_ind, x_target, x_final<jupyter_output><empty_output><jupyter_text>Collect experiment results using the read_finalresult function- collect final error (using ADG and MSE functions defined before)- collect initial c and x to computer initial error (using ADG and MSE functions defined before)<jupyter_code>def Collect_Err(num_lp,  nVar, num_constraints,loss):
    err_final = np.zeros((num_lp,len(num_constraints)))
    err_init = np.zeros((num_lp, len(num_constraints)))
    for j,nCons in enumerate (num_constraints): 
        print("LPs with %dVar, %dCons"%(nVar,nCons))
        result, c_init, c_final, best_ind, x_target, x_final = read_finalresult(num_lp, nCons, nVar, loss)
        err = result[best_ind,-2]
        err_final[:,j] = np.array(err).ravel()
        
        for ind_lp in range(num_lp):
            filename = PATH + "Deep_Inv_Opt/figure3_i/n=%d,m=%d/%svar_LP%s.txt"%(nVar, nCons, str(nVar).zfill(2),str(ind_lp+1).zfill(2))
            A, b, _, _ = Read_3aLP(filename)
            
            A_ub = io.tensor(A)
            b_ub = io.tensor(b).view(-1,1)
            c_vector = io.tensor(c_init[ind_lp]).view(-1,1)
            c_vector/=sum(abs(c_vector))
            x_predict = io.linprog(c_vector, A_ub, b_ub, eps = 1e-6).detach().numpy()
            if loss == "ADG":            
                err_init[ind_lp, j] = ADG( c_vector.detach().numpy(), 
                                        x_predict, x_target[ind_lp])
            elif loss == "MSE":            
                err_init[ind_lp, j] = MSE(x_predict, x_target[ind_lp])
    return err_init, err_final

num_lp = 50
RANGE_CONS = [[4,  8,  16],
               [20, 36, 80]]

ADG02_init, ADG02_final = Collect_Err(num_lp, 2, RANGE_CONS[0], "ADG")

ADG10_init, ADG10_final = Collect_Err(num_lp, 10, RANGE_CONS[1], "ADG")

MSE02_init, MSE02_final = Collect_Err(num_lp, 2, RANGE_CONS[0], "MSE")

MSE10_init, MSE10_final = Collect_Err(num_lp, 10, RANGE_CONS[1], "MSE")<jupyter_output><empty_output><jupyter_text>Create boxplot using the collected data from the previous cell**Note, the final boxplots might not look identical to what we presented in the paper due to the randomization in baseline LP generation.**<jupyter_code>MSE_ERR = (MSE02_init, MSE02_final, MSE10_init, MSE10_final)
ADG_ERR = (ADG02_init, ADG02_final, ADG10_init, ADG10_final )

def final_plot(err_data, loss):
    init02, final02, init10, final10 = err_data
    print(final02.shape)
        
    data02 = init02
    for i in range(3):
        data02 = np.insert(data02, i*2+1, final02[:,i], axis=1)
    data10 = init10
    for i in range(3):
        data10 = np.insert(data10, i*2+1, final10[:,i], axis=1)
        
        
    data_ = np.column_stack((data02,data10))

    if _LOSS == "ADG":
        data_ [data_<1e-8]=1e-8
    elif _LOSS == "MSE":
        data_ [data_<1e-10]=1e-10
    
    return data_

for _LOSS in ["MSE", "ADG"]:
    if _LOSS == "ADG":
        data = final_plot(ADG_ERR, _LOSS)
        fig = plt.figure(figsize=(10,10))
        ax = fig.add_subplot(111)
        plt.ylim((1e-9,1e1))
    elif _LOSS == "MSE":
        data = final_plot(MSE_ERR, _LOSS)
        fig = plt.figure(figsize=(10,10))
        ax = fig.add_subplot(111)
        plt.ylim((1e-11,1e2))

    _, nbox = data.shape    

    ind_init = np.arange(0,nbox,2)
    ind_final = np.arange(1,nbox,2)

    y_init = data[:, ind_init]
    x_init = np.random.uniform(1-.15, 1+0.15, size=len(data[:, ind_init]))
    y_final = data[:, ind_final]
    x_final = np.random.uniform(2-.15, 2+0.15, size=len(data[:, ind_final]))
    for i in np.arange(2, nbox,2):
        min = i+1-0.15
        max = i+1+0.15
        x_init = np.column_stack((x_init, np.random.uniform(i+1-.15, i+1+0.15, 
                                                            size=len(data[:, i]))))
    for i in np.arange(3, nbox,2):
        min = i+1-0.15
        max = i+1+0.15
        x_final = np.column_stack((x_final, np.random.uniform(i+1-.15, i+1+0.15, 
                                                              size=len(data[:, i]))))
    ax.boxplot(data,sym='',medianprops =dict(color='black'))
    ax.plot([], 'go',markerfacecolor='w', alpha=0.8, label = 'init_loss vertex target')  
    ax.plot([], 'bo',markerfacecolor='w', alpha=0.8, label = 'final_loss vertex target') 

    for i in range(nbox//2):
        ax.plot(x_init[:, i], y_init[:, i], 'go',markerfacecolor='w', alpha=0.8) 

    for i in range(nbox//2):
        ax.plot(x_final[:, i], y_final[:, i], 'bo',markerfacecolor='w', alpha=0.8) 
    plt.legend(loc='upper center', bbox_to_anchor=(0.5, -0.1),
              fontsize = 12, ncol=3)

    
    plt.xticks([], [])
    plt.yscale("log")
    plt.yticks(fontsize=14)
    plt.savefig(PATH+"Deep_Inv_Opt/deepinverse_figure3_i_%s.pdf"%_LOSS, frameon = True, bbox_inches='tight', dpi=100)
    plt.show()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14289237446976422, 'fraction_numerical': 0.02337969813554306, 'mean_word_length': 2.408135821146411, 'pattern_counts': {'":': 10, '<': 42, '<?xml version=': 0, '>': 40, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 30, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '9', 'dir': 'jupyter-structured-clean-dedup', 'id': '419965', 'license': 'MIT', 'n_tokens_mistral': 7378, 'n_tokens_neox': 7088, 'n_words': 1602, 'path': 'paper/Figure 3_step2_(i)_Learn c.ipynb', 'repo_name': 'tankconcordia/deep_inv_opt', 'size': '20295'}