Warmstarts.

.gitignore

@@ -24,4 +24,5 @@ __pycache__
 *.tree
 
 
-/dump/
+/dump/
+tatus

data/baron/minlp.jl

@@ -11,7 +11,7 @@
 # // Objective nonlinear
 # // Nonlinear constraints
 
-function minlp(gm::Bool = false, eps = 1e-5)
+function minlp(gm::Bool = false, eps = 1e-8)
     m = JuMP.Model()
     @variable(m, 0 <= y[1:3] <= 1, Int)
     @variable(m, x[1:3] >= eps)

src/OptimalConstraintTree.jl

@@ -12,8 +12,9 @@ module OptimalConstraintTree
 #     using Plots
     using ProgressMeter
     using Random
+    using StatsBase
 
-    const PROJECT_ROOT = pwd()
+    const PROJECT_ROOT = dirname(dirname(@__FILE__))
     const SAVE_DIR = PROJECT_ROOT * "\\data\\"
     const GAMS_DIR = PROJECT_ROOT * "\\data\\gams\\"
 

src/fit.jl

@@ -21,7 +21,7 @@ function lnr_kwargs(; kwargs...)
     # TODO: figure out how to merge this stuff!!!!
     nkwargs = Dict{Symbol, Any}() 
     valid_keys = [:random_seed, :max_depth, :cp, :minbucket, :fast_num_support_restarts, 
-                  :localsearch, :ls_num_hyper_restarts, :ls_num_tree_restarts]
+                  :localsearch, :ls_num_hyper_restarts, :ls_num_tree_restarts, :hyperplane_config]
     for item in kwargs
         if item.first in valid_keys
             nkwargs[item.first] = item.second

src/global_model.jl

@@ -367,7 +367,11 @@ function accuracy(bbf::Union{GlobalModel, BlackBoxFunction})
     end
 end
 
-""" Extends JuMP.optimize! to GlobalModels. """
+""" 
+    JuMP.optimize!(gm::GlobalModel; kwargs...)
+
+Applies JuMP.optimize! to GlobalModels, and saves solution history. 
+"""
 function JuMP.optimize!(gm::GlobalModel; kwargs...)
     JuMP.optimize!(gm.model, kwargs...)
     append!(gm.solution_history, solution(gm), cols=:setequal)
@@ -433,6 +437,8 @@ end
 
 match_bbfs_to_vars(gm::GlobalModel, vars::Array = JuMP.all_variables(gm)) = match_bbfs_to_vars(gm.bbfs, vars)
 
+""" Clears all sampling, training and optimization data from GlobalModel."""
 function clear_data!(gm::GlobalModel)
     clear_data!.(gm.bbfs)
+    gm.solution_history = DataFrame([Float64 for i=1:length(gm.vars)], string.(gm.vars))
 end

src/learners.jl

@@ -17,12 +17,11 @@ function base_otr()
 end
 
 """
-Returns the baseline OptimalTreeClassifier.
+Returns the baseline OptimalTreeClassifier, with all modifyable args. 
 """
 function base_otc()
     lnr = IAI.OptimalTreeClassifier(
         hyperplane_config = (sparsity = :all,),
-        # Modifiables, with defaults
         random_seed = 1,
         max_depth = 5, 
         cp = 1e-6,

src/plot.jl

@@ -1,16 +1,3 @@
-function plot(tree::IAI.GridSearch)
-    """ Shows grid in browser. """
-    try
-        IAI.show_in_browser(tree.lnr)
-    catch err
-        if isa(err, UndefRefError)
-        @warn("Certain trees are untrained.")
-        else
-            rethrow(err)
-        end
-    end
-end
-
 function show_trees(bbf::Union{BlackBoxFunction, DataConstraint})
     for grid in bbf.learners
         IAI.show_in_browser(grid.lnr)

src/warmstart.jl

@@ -0,0 +1,45 @@
+""" Warmstarts an OptimalTreeLearner using previous solutions. """
+function fit_warmstart!(lnr, X, y, tree)
+  IAIBase.resetinternalparams!(lnr)
+  criterion = IAIBase.process_criterion(lnr)
+  IAIBase.validate_criterion(criterion, lnr, :criterion)
+  data = IAIBase.Data(tree, X, nothing, y)
+  mins, maxs = IAIBase.get_extrema(data.features.X_numeric,
+                                   data.features.n_numeric_orig_features)
+  data.features.X_mins .= mins
+  data.features.X_maxs .= maxs
+  lnr.normalize_X && IAIBase.normalize_X!(data.features)
+  IAIBase._preprocess!(lnr, data, IAIBase.gettasktype(lnr))
+  lnr.prb_ = IAIBase.Problem(lnr, data, criterion)
+
+  ls = OptimalTrees.LocalSearcher(lnr, lnr.prb_)
+  OptimalTrees.init!(ls, lnr)
+  Random.seed!(ls.rng_gen, ls.random_seed)
+
+  lnr.all_trees_ = deepcopy(tree.all_trees_)
+  for t in 1:length(lnr.all_trees_)
+    rng = IAIBase.RandomStreams.next_stream(ls.rng_gen)
+    OptimalTrees.set_rng!(ls, rng)
+
+    tree = lnr.all_trees_[t]
+
+    if lnr.normalize_X
+      IAIBase.normalize_X!(tree, lnr.prb_.data.features)
+    end
+    if :normalize_y in fieldnames(typeof(lnr)) && lnr.normalize_y
+      IAIBase.normalize_y!(tree, lnr.prb_.data.target)
+    end
+
+    OptimalTrees.local_search!(tree, ls)
+
+    IAITrees.finalize_capacity!(tree)
+    OptimalTrees.update_tree_fits!(tree, ls.data, ls.evaluators, true)
+    lnr.all_trees_[t] = tree
+
+  end
+  best_ind = argmin([OptimalTrees.train_error(lnr, tree)
+                     for tree in lnr.all_trees_])
+  lnr.tree_ = lnr.all_trees_[best_ind]
+
+  lnr
+end

test/afpm.jl

@@ -3,16 +3,15 @@ include("../test/load.jl")
 using DataFrames
 using Test
 using CSV
-using Gurobi
 
 @info "Trains over axial flux motor data."
 
 # Getting simulation data
-X = CSV.read("data/afpm/afpm_inputs.csv",
+X = CSV.read("data/afpm/afpm_inputs.csv", DataFrame,
                  copycols=true, delim=",");
-X_infeas = CSV.read("data/afpm/afpm_infeas_inputs.csv",
+X_infeas = CSV.read("data/afpm/afpm_infeas_inputs.csv", DataFrame,
                     copycols=true, delim=",");
-Y = CSV.read("data/afpm/afpm_outputs.csv",
+Y = CSV.read("data/afpm/afpm_outputs.csv", DataFrame,
                  copycols=true, delim=",");
 X = select!(X, Not(:Column1));
 X_infeas = select!(X_infeas, Not(:Column1));
@@ -23,7 +22,7 @@ varkeys = Symbol.(names(X))
 objkeys = Symbol.(names(Y))
 
 # Creating model with relevant variables
-m = Model(Gurobi.Optimizer)
+m = Model(with_optimizer(CPLEX_SILENT))
 inputs = [];
 for key in varkeys
     nv = @variable(m, base_name = string(key))
@@ -43,33 +42,33 @@ N_coils, TPC, p = inputs[4:6];
 # ranges = Dict(key => [0.5*val, 1.5*val] for (key, val) in baseline);
 # ranges[inputs[9]] = [0.25*baseline[inputs[9]], 2.25*baseline[inputs[9]]]; # wire_A needs special care...
 # ranges = Dict(key => log.(val) for (key, val) in ranges)
-ranges = Dict(var => [log.(minimum(X[string(var)])),log.(maximum(X[string(var)]))] for var in inputs)
+ranges = Dict(var => [log.(minimum(X[!, string(var)])),log.(maximum(X[!, string(var)]))] for var in inputs)
 
 # Bounding for output variables of interest
 idxs = [1, 4, 6, 8, 10, 12];
 feas_idxs = findall(x -> x .>= 0.5, Y.Efficiency);
 X_feas = X[feas_idxs, :];
 Y_feas = Y[feas_idxs, :];
-out_ranges = Dict(key => [minimum((Y_feas[Symbol(key)])), maximum((Y_feas[Symbol(key)]))] for key in outputs[idxs]);
+out_ranges = Dict(key => [minimum((Y_feas[!, Symbol(key)])), maximum((Y_feas[!, Symbol(key)]))] for key in outputs[idxs]);
 bound!(m, ranges);
 
 # Geometry constraints (in logspace)
 D_out, D_in, N_coils, wire_w = inputs[1], inputs[2], inputs[4], inputs[7]
 @constraint(m, D_out >= D_in)
 @constraint(m, log(pi) + D_in >= log(0.2) + wire_w + N_coils) #pi*D_in >= 2*0.1*wire_w*N_coils
 
-N_coils_range = log.(unique(X_feas["N_coils"]));
+N_coils_range = log.(unique(X_feas[!, "N_coils"]));
 int = @variable(m, [1:length(N_coils_range)], Bin)
 @constraint(m, sum(int) == 1)
 @constraint(m, N_coils == sum(N_coils_range .* int))
 
-p_range = log.(unique(X_feas["p"]))
+p_range = log.(unique(X_feas[!, "p"]))
 int = @variable(m, [1:length(p_range)], Bin)
 @constraint(m, sum(int) == 1)
 @constraint(m, p == sum(p_range .* int))
 @constraint(m, N_coils >= p + 1e-3) # motor type 2
 
-TPC_range = log.(unique(X_feas["TPC"]))
+TPC_range = log.(unique(X_feas[!, "TPC"]))
 int = @variable(m, [1:length(TPC_range)], Bin)
 @constraint(m, sum(int) == 1)
 @constraint(m, TPC == sum(TPC_range .* int))
@@ -80,8 +79,10 @@ P_shaft, Torque, Rotational_Speed, Efficiency, Mass, Mass_Specific_Power = [outp
 set_upper_bound(Efficiency, 1)
 
 # Fitting power closure, and creating a global model
-# simulation = DataConstraint(vars = inputs)
 feasmap = zeros(size(Y, 1)); feasmap[feas_idxs] .= 1;
+# simulation = DataConstraint(xvars = inputs, yvars = nothing)
+# simulation.X = log.(X)
+# simulation.Y = feasmap
 # add_data!(simulation, log.(X), feasmap)
 # learn_constraint!(simulation)
 # lnr = IAI.fit!(base_otc(), log.(X), feasmap)
@@ -97,7 +98,7 @@ P_shaft, Torque, Rotational_Speed, Efficiency, Mass, Mass_Specific_Power = [outp
 FOMs = [P_shaft, Rotational_Speed, Mass, Efficiency]
 leaf_index, all_leaves = bin_to_leaves(lnr, log.(X_feas))
 for FOM in FOMs
-    regressor = regress(log.(X_feas), log.(Y_feas[string(FOM)]))
+    regressor = regress(log.(X_feas), log.(Y_feas[!,string(FOM)]))
     constant = IAI.get_prediction_constant(regressor)
     weights  = IAI.get_prediction_weights(regressor)[1]
     vks = Symbol.(names(X_feas))
@@ -117,14 +118,14 @@ end
 @constraint(m, log.(7800) <= Rotational_Speed <= log(8200))
 @objective(m, Min, Mass)
 optimize!(m)
-# println("Inputs")
-# for i=1:length(inputs)
-#     println(string(inputs[i], " ", exp(getvalue(inputs[i]))))
-# end
-# println("FOMs")
-# for FOM in FOMs
-#     println(string(FOM, " ", exp(getvalue(FOM))))
-# end
+println("Inputs")
+for i=1:length(inputs)
+    println(string(inputs[i], " ", exp(getvalue(inputs[i]))))
+end
+println("FOMs")
+for FOM in FOMs
+    println(string(FOM, " ", exp(getvalue(FOM))))
+end
 
 # fdf = DataFrame(names(X) .=> exp.(getvalue.(inputs)))
 # CSV.write("data/afpm/afpm_opt.csv", fdf)

test/all.jl

@@ -3,9 +3,10 @@ all_tests:
 - Julia version: 
 - Author: Berk
 - Date: 2020-07-01
-Note: All tests should be run from a julia REPR within the julia folder, using:
+Note: All tests should be run from a julia REPR within the OptimalConstraintTree folder, using:
       julia --project=.
       include("test/all.jl")
+
       To see coverage, run with:
       julia --project=. --code-coverage=tracefile-%p.info --code-coverage=user
       include("test/all.jl")