ANL-CEEESA/MIPLearn
1
0
Fork 0
mirror of https://github.com/ANL-CEEESA/MIPLearn.git synced 2025-12-06 09:28:51 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'master' into feature/tsp

Browse Source
This commit is contained in:
Alinson S. Xavier
2020-02-24 14:22:20 -06:00
parent 2dc66696f9 f22aecb14d
commit f713a399a8
107 changed files with 13071 additions and 914 deletions
Download Patch File Download Diff File Expand all files Collapse all files

25
miplearn/__init__.py
View File

@@ -1,12 +1,19 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from .component import Component
from .extractors import (SolutionExtractor,
CombinedExtractor,
InstanceFeaturesExtractor,
ObjectiveValueExtractor,
VariableFeaturesExtractor,
)
from .components.component import Component
from .components.objective import ObjectiveValueComponent
from .components.primal import (PrimalSolutionComponent,
AdaptivePredictor,
)
from .components.branching import BranchPriorityComponent
from .benchmark import BenchmarkRunner
from .instance import Instance
from .solvers import LearningSolver
from .benchmark import BenchmarkRunner
from .warmstart import WarmStartComponent, KnnWarmStartPredictor, LogisticWarmStartPredictor
from .branching import BranchPriorityComponent
from .extractors import UserFeaturesExtractor, SolutionExtractor

109
miplearn/benchmark.py
View File

@@ -1,9 +1,11 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from .solvers import LearningSolver
from copy import deepcopy
import pandas as pd
from tqdm.auto import tqdm
class BenchmarkRunner:
def __init__(self, solvers):
@@ -13,59 +15,23 @@ class BenchmarkRunner:
self.solvers = solvers
self.results = None
def solve(self, instances, fit=True, tee=False):
for (name, solver) in self.solvers.items():
for i in tqdm(range(len((instances)))):
results = solver.solve(deepcopy(instances[i]), tee=tee)
self._push_result(results, solver=solver, name=name, instance=i)
if fit:
solver.fit()
def parallel_solve(self, instances, n_jobs=1, n_trials=1):
if self.results is None:
self.results = pd.DataFrame(columns=["Solver",
"Instance",
"Wallclock Time",
"Lower Bound",
"Upper Bound",
"Gap",
"Nodes",
])
instances = instances * n_trials
for (name, solver) in self.solvers.items():
results = solver.parallel_solve(instances,
n_jobs=n_jobs,
label=name,
label="Solve (%s)" % name,
collect_training_data=False)
for i in range(len(instances)):
wallclock_time = None
for key in ["Time", "Wall time", "Wallclock time"]:
if key not in results[i]["Solver"][0].keys():
continue
if str(results[i]["Solver"][0][key]) == "<undefined>":
continue
wallclock_time = float(results[i]["Solver"][0][key])
nodes = results[i]["Solver"][0]["Nodes"]
lb = results[i]["Problem"][0]["Lower bound"]
ub = results[i]["Problem"][0]["Upper bound"]
gap = (ub - lb) / lb
self.results = self.results.append({
"Solver": name,
"Instance": i,
"Wallclock Time": wallclock_time,
"Lower Bound": lb,
"Upper Bound": ub,
"Gap": gap,
"Nodes": nodes,
}, ignore_index=True)
groups = self.results.groupby("Instance")
best_lower_bound = groups["Lower Bound"].transform("max")
best_upper_bound = groups["Upper Bound"].transform("min")
best_gap = groups["Gap"].transform("min")
best_nodes = groups["Nodes"].transform("min")
best_wallclock_time = groups["Wallclock Time"].transform("min")
self.results["Relative Lower Bound"] = \
self.results["Lower Bound"] / best_lower_bound
self.results["Relative Upper Bound"] = \
self.results["Upper Bound"] / best_upper_bound
self.results["Relative Wallclock Time"] = \
self.results["Wallclock Time"] / best_wallclock_time
self.results["Relative Gap"] = \
self.results["Gap"] / best_gap
self.results["Relative Nodes"] = \
self.results["Nodes"] / best_nodes
self._push_result(results[i], solver=solver, name=name, instance=i)
def raw_results(self):
return self.results
@@ -80,6 +46,47 @@ class BenchmarkRunner:
for (name, solver) in self.solvers.items():
solver.load_state(filename)
def fit(self):
def fit(self, training_instances):
for (name, solver) in self.solvers.items():
solver.fit()
solver.fit(training_instances)
def _push_result(self, result, solver, name, instance):
if self.results is None:
self.results = pd.DataFrame(columns=["Solver",
"Instance",
"Wallclock Time",
"Lower Bound",
"Upper Bound",
"Gap",
"Nodes",
"Mode",
])
lb = result["Lower bound"]
ub = result["Upper bound"]
gap = (ub - lb) / lb
self.results = self.results.append({
"Solver": name,
"Instance": instance,
"Wallclock Time": result["Wallclock time"],
"Lower Bound": lb,
"Upper Bound": ub,
"Gap": gap,
"Nodes": result["Nodes"],
"Mode": solver.mode,
}, ignore_index=True)
groups = self.results.groupby("Instance")
best_lower_bound = groups["Lower Bound"].transform("max")
best_upper_bound = groups["Upper Bound"].transform("min")
best_gap = groups["Gap"].transform("min")
best_nodes = groups["Nodes"].transform("min")
best_wallclock_time = groups["Wallclock Time"].transform("min")
self.results["Relative Lower Bound"] = \
self.results["Lower Bound"] / best_lower_bound
self.results["Relative Upper Bound"] = \
self.results["Upper Bound"] / best_upper_bound
self.results["Relative Wallclock Time"] = \
self.results["Wallclock Time"] / best_wallclock_time
self.results["Relative Gap"] = \
self.results["Gap"] / best_gap
self.results["Relative Nodes"] = \
self.results["Nodes"] / best_nodes

23
miplearn/component.py
View File

@@ -1,23 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
from abc import ABC, abstractmethod
class Component(ABC):
@abstractmethod
def fit(self, solver):
pass
@abstractmethod
def before_solve(self, solver, instance, model):
pass
@abstractmethod
def after_solve(self, solver, instance, model):
pass
@abstractmethod
def merge(self, other):
pass

4
miplearn/components/__init__.py Normal file
View File

@@ -0,0 +1,4 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

4
miplearn/scripts/branchpriority.jl → miplearn/components/branching.jl
View File

@@ -1,3 +1,7 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
import Base.Threads.@threads
using TinyBnB, CPLEXW, Printf

34
miplearn/branching.py → miplearn/components/branching.py
View File

@@ -1,34 +1,37 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from . import Component
from .transformers import PerVariableTransformer
from .component import Component
from ..extractors import Extractor
from abc import ABC, abstractmethod
from sklearn.neighbors import KNeighborsRegressor
import numpy as np
from p_tqdm import p_map
from tqdm.auto import tqdm
from joblib import Parallel, delayed
import multiprocessing
def _default_branch_priority_predictor():
return KNeighborsRegressor(n_neighbors=1)
class BranchPriorityComponent(Component):
def __init__(self,
node_limit=1_000,
node_limit=10_000,
predictor=_default_branch_priority_predictor,
):
self.transformer = PerVariableTransformer()
self.pending_instances = []
self.x_train = {}
self.y_train = {}
self.predictors = {}
self.node_limit = node_limit
self.predictor_factory = predictor
def before_solve(self, solver, instance, model):
assert solver.is_persistent, "BranchPriorityComponent requires a persistent solver"
assert solver.internal_solver.name == "gurobi_persistent", "Only GurobiPersistent is currently supported"
from gurobipy import GRB
var_split = self.transformer.split_variables(instance, model)
var_split = Extractor.split_variables(instance, model)
for category in var_split.keys():
if category not in self.predictors.keys():
continue
@@ -56,7 +59,7 @@ class BranchPriorityComponent(Component):
src_dirname = os.path.dirname(os.path.realpath(__file__))
priority_file = tempfile.NamedTemporaryFile(mode="r")
subprocess.run(["julia",
"%s/scripts/branchpriority.jl" % src_dirname,
"%s/branching.jl" % src_dirname,
lp_file.name,
priority_file.name,
str(self.node_limit),
@@ -96,7 +99,7 @@ class BranchPriorityComponent(Component):
return comp
# Run strong branching on pending instances
subcomponents = Parallel(n_jobs=n_jobs)(
delayed(_process)(instance)
for instance in tqdm(self.pending_instances, desc="Branch priority")
@@ -108,9 +111,8 @@ class BranchPriorityComponent(Component):
for category in self.x_train.keys():
x_train = self.x_train[category]
y_train = self.y_train[category]
self.predictors[category] = KNeighborsRegressor(n_neighbors=1)
self.predictors[category] = self.predictor_factory()
self.predictors[category].fit(x_train, y_train)
def _build_x(self, instance, var, index):
instance_features = instance.get_instance_features()

27
miplearn/components/component.py Normal file
View File

@@ -0,0 +1,27 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from abc import ABC, abstractmethod
class Component(ABC):
"""
A Component is an object which adds functionality to a LearningSolver.
"""
@abstractmethod
def before_solve(self, solver, instance, model):
pass
@abstractmethod
def after_solve(self, solver, instance, model):
pass
@abstractmethod
def merge(self, other):
pass
@abstractmethod
def fit(self, training_instances):
pass

49
miplearn/components/objective.py Normal file
View File

@@ -0,0 +1,49 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from .. import Component, InstanceFeaturesExtractor, ObjectiveValueExtractor
from sklearn.linear_model import LinearRegression
from copy import deepcopy
import numpy as np
import logging
logger = logging.getLogger(__name__)
class ObjectiveValueComponent(Component):
"""
A Component which predicts the optimal objective value of the problem.
"""
def __init__(self,
regressor=LinearRegression()):
self.ub_regressor = None
self.lb_regressor = None
self.regressor_prototype = regressor
def before_solve(self, solver, instance, model):
if self.ub_regressor is not None:
lb, ub = self.predict([instance])[0]
instance.predicted_ub = ub
instance.predicted_lb = lb
logger.info("Predicted objective: [%.2f, %.2f]" % (lb, ub))
def after_solve(self, solver, instance, model):
pass
def merge(self, other):
pass
def fit(self, training_instances):
features = InstanceFeaturesExtractor().extract(training_instances)
ub = ObjectiveValueExtractor(kind="upper bound").extract(training_instances)
lb = ObjectiveValueExtractor(kind="lower bound").extract(training_instances)
self.ub_regressor = deepcopy(self.regressor_prototype)
self.lb_regressor = deepcopy(self.regressor_prototype)
self.ub_regressor.fit(features, ub)
self.lb_regressor.fit(features, lb)
def predict(self, instances):
features = InstanceFeaturesExtractor().extract(instances)
lb = self.lb_regressor.predict(features)
ub = self.ub_regressor.predict(features)
return np.hstack([lb, ub])

205
miplearn/components/primal.py Normal file
View File

@@ -0,0 +1,205 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from .component import Component
from ..extractors import *
from abc import ABC, abstractmethod
from copy import deepcopy
import numpy as np
from sklearn.pipeline import make_pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import cross_val_score
from sklearn.metrics import roc_curve
from sklearn.neighbors import KNeighborsClassifier
from tqdm.auto import tqdm
import pyomo.environ as pe
import logging
logger = logging.getLogger(__name__)
class AdaptivePredictor:
def __init__(self,
predictor=None,
min_samples_predict=1,
min_samples_cv=100,
thr_fix=0.999,
thr_alpha=0.50,
thr_balance=0.95,
):
self.min_samples_predict = min_samples_predict
self.min_samples_cv = min_samples_cv
self.thr_fix = thr_fix
self.thr_alpha = thr_alpha
self.thr_balance = thr_balance
self.predictor_factory = predictor
def fit(self, x_train, y_train):
n_samples = x_train.shape[0]
# If number of samples is too small, don't predict anything.
if n_samples < self.min_samples_predict:
logger.debug(" Too few samples (%d); always predicting false" % n_samples)
self.predictor = 0
return
# If vast majority of observations are false, always return false.
y_train_avg = np.average(y_train)
if y_train_avg <= 1.0 - self.thr_fix:
logger.debug(" Most samples are negative (%.3f); always returning false" % y_train_avg)
self.predictor = 0
return
# If vast majority of observations are true, always return true.
if y_train_avg >= self.thr_fix:
logger.debug(" Most samples are positive (%.3f); always returning true" % y_train_avg)
self.predictor = 1
return
# If classes are too unbalanced, don't predict anything.
if y_train_avg < (1 - self.thr_balance) or y_train_avg > self.thr_balance:
logger.debug(" Classes are too unbalanced (%.3f); always returning false" % y_train_avg)
self.predictor = 0
return
# Select ML model if none is provided
if self.predictor_factory is None:
if n_samples < 30:
self.predictor_factory = KNeighborsClassifier(n_neighbors=n_samples)
else:
self.predictor_factory = make_pipeline(StandardScaler(), LogisticRegression())
# Create predictor
if callable(self.predictor_factory):
pred = self.predictor_factory()
else:
pred = deepcopy(self.predictor_factory)
# Skip cross-validation if number of samples is too small
if n_samples < self.min_samples_cv:
logger.debug(" Too few samples (%d); skipping cross validation" % n_samples)
self.predictor = pred
self.predictor.fit(x_train, y_train)
return
# Calculate cross-validation score
cv_score = np.mean(cross_val_score(pred, x_train, y_train, cv=5))
dummy_score = max(y_train_avg, 1 - y_train_avg)
cv_thr = 1. * self.thr_alpha + dummy_score * (1 - self.thr_alpha)
# If cross-validation score is too low, don't predict anything.
if cv_score < cv_thr:
logger.debug(" Score is too low (%.3f < %.3f); always returning false" % (cv_score, cv_thr))
self.predictor = 0
else:
logger.debug(" Score is acceptable (%.3f > %.3f); training classifier" % (cv_score, cv_thr))
self.predictor = pred
self.predictor.fit(x_train, y_train)
def predict_proba(self, x_test):
if isinstance(self.predictor, int):
y_pred = np.zeros((x_test.shape[0], 2))
y_pred[:, self.predictor] = 1.0
return y_pred
else:
return self.predictor.predict_proba(x_test)
class PrimalSolutionComponent(Component):
"""
A component that predicts primal solutions.
"""
def __init__(self,
predictor=AdaptivePredictor(),
mode="exact",
max_fpr=[1e-3, 1e-3],
min_threshold=[0.75, 0.75],
dynamic_thresholds=True,
):
self.mode = mode
self.predictors = {}
self.is_warm_start_available = False
self.max_fpr = max_fpr
self.min_threshold = min_threshold
self.thresholds = {}
self.predictor_factory = predictor
self.dynamic_thresholds = dynamic_thresholds
def before_solve(self, solver, instance, model):
solution = self.predict(instance, model)
if self.mode == "heuristic":
solver.internal_solver.fix(solution)
else:
solver.internal_solver.set_warm_start(solution)
def after_solve(self, solver, instance, model):
pass
def fit(self, training_instances):
features = VariableFeaturesExtractor().extract(training_instances)
solutions = SolutionExtractor().extract(training_instances)
for category in tqdm(features.keys(), desc="Fit (Primal)"):
x_train = features[category]
y_train = solutions[category]
for label in [0, 1]:
logger.debug("Fitting predictors[%s, %s]:" % (category, label))
if callable(self.predictor_factory):
pred = self.predictor_factory(category, label)
else:
pred = deepcopy(self.predictor_factory)
self.predictors[category, label] = pred
y = y_train[:, label].astype(int)
pred.fit(x_train, y)
# If y is either always one or always zero, set fixed threshold
y_avg = np.average(y)
if (not self.dynamic_thresholds) or y_avg <= 0.001 or y_avg >= 0.999:
self.thresholds[category, label] = self.min_threshold[label]
logger.debug(" Setting threshold to %.4f" % self.min_threshold[label])
continue
# Calculate threshold dynamically using ROC curve
y_scores = pred.predict_proba(x_train)[:, 1]
fpr, tpr, thresholds = roc_curve(y, y_scores)
k = 0
while True:
if (k + 1) > len(fpr):
break
if fpr[k + 1] > self.max_fpr[label]:
break
if thresholds[k + 1] < self.min_threshold[label]:
break
k = k + 1
logger.debug(" Setting threshold to %.4f (fpr=%.4f, tpr=%.4f)"%
(thresholds[k], fpr[k], tpr[k]))
self.thresholds[category, label] = thresholds[k]
def predict(self, instance, model=None):
if model is None:
model = instance.to_model()
x_test = VariableFeaturesExtractor().extract([instance], [model])
solution = {}
var_split = Extractor.split_variables(instance, model)
for category in var_split.keys():
for (i, (var, index)) in enumerate(var_split[category]):
if var not in solution.keys():
solution[var] = {}
solution[var][index] = None
for label in [0, 1]:
if (category, label) not in self.predictors.keys():
continue
ws = self.predictors[category, label].predict_proba(x_test[category])
logger.debug("%s[%s] ws=%.6f threshold=%.6f" %
(var, index, ws[i, 1], self.thresholds[category, label]))
if ws[i, 1] >= self.thresholds[category, label]:
solution[var][index] = label
return solution
def merge(self, other_components):
pass

4
miplearn/components/tests/__init__.py Normal file
View File

@@ -0,0 +1,4 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

49
miplearn/components/tests/test_branching.py Normal file
View File

@@ -0,0 +1,49 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import BranchPriorityComponent, LearningSolver
from miplearn.problems.knapsack import KnapsackInstance
import numpy as np
import tempfile
def _get_instances():
return [
KnapsackInstance(
weights=[23., 26., 20., 18.],
prices=[505., 352., 458., 220.],
capacity=67.,
),
] * 2
def test_branching():
instances = _get_instances()
component = BranchPriorityComponent()
for instance in instances:
component.after_solve(None, instance, None)
component.fit(None)
for key in ["default"]:
assert key in component.x_train.keys()
assert key in component.y_train.keys()
assert component.x_train[key].shape == (8, 4)
assert component.y_train[key].shape == (8, 1)
# def test_branch_priority_save_load():
# state_file = tempfile.NamedTemporaryFile(mode="r")
# solver = LearningSolver(components={"branch-priority": BranchPriorityComponent()})
# solver.parallel_solve(_get_instances(), n_jobs=2)
# solver.fit()
# comp = solver.components["branch-priority"]
# assert comp.x_train["default"].shape == (8, 4)
# assert comp.y_train["default"].shape == (8, 1)
# assert "default" in comp.predictors.keys()
# solver.save_state(state_file.name)
#
# solver = LearningSolver(components={"branch-priority": BranchPriorityComponent()})
# solver.load_state(state_file.name)
# comp = solver.components["branch-priority"]
# assert comp.x_train["default"].shape == (8, 4)
# assert comp.y_train["default"].shape == (8, 1)
# assert "default" in comp.predictors.keys()

29
miplearn/components/tests/test_objective.py Normal file
View File

@@ -0,0 +1,29 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import ObjectiveValueComponent, LearningSolver
from miplearn.problems.knapsack import KnapsackInstance
def _get_instances():
instances = [
KnapsackInstance(
weights=[23., 26., 20., 18.],
prices=[505., 352., 458., 220.],
capacity=67.,
),
]
models = [instance.to_model() for instance in instances]
solver = LearningSolver()
for i in range(len(instances)):
solver.solve(instances[i], models[i])
return instances, models
def test_usage():
instances, models = _get_instances()
comp = ObjectiveValueComponent()
comp.fit(instances)
assert instances[0].lower_bound == 1183.0
assert instances[0].upper_bound == 1183.0
assert comp.predict(instances).tolist() == [[1183.0, 1183.0]]

55
miplearn/components/tests/test_primal.py Normal file
View File

@@ -0,0 +1,55 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import LearningSolver, PrimalSolutionComponent
from miplearn.problems.knapsack import KnapsackInstance
import numpy as np
import tempfile
def _get_instances():
instances = [
KnapsackInstance(
weights=[23., 26., 20., 18.],
prices=[505., 352., 458., 220.],
capacity=67.,
),
] * 5
models = [inst.to_model() for inst in instances]
solver = LearningSolver()
for i in range(len(instances)):
solver.solve(instances[i], models[i])
return instances, models
def test_predict():
instances, models = _get_instances()
comp = PrimalSolutionComponent()
comp.fit(instances)
solution = comp.predict(instances[0], models[0])
assert models[0].x in solution.keys()
for idx in range(4):
assert idx in solution[models[0].x].keys()
# def test_warm_start_save_load():
# state_file = tempfile.NamedTemporaryFile(mode="r")
# solver = LearningSolver(components={"warm-start": WarmStartComponent()})
# solver.parallel_solve(_get_instances(), n_jobs=2)
# solver.fit()
# comp = solver.components["warm-start"]
# assert comp.x_train["default"].shape == (8, 6)
# assert comp.y_train["default"].shape == (8, 2)
# assert ("default", 0) in comp.predictors.keys()
# assert ("default", 1) in comp.predictors.keys()
# solver.save_state(state_file.name)
# solver.solve(_get_instances()[0])
# solver = LearningSolver(components={"warm-start": WarmStartComponent()})
# solver.load_state(state_file.name)
# comp = solver.components["warm-start"]
# assert comp.x_train["default"].shape == (8, 6)
# assert comp.y_train["default"].shape == (8, 2)
# assert ("default", 0) in comp.predictors.keys()
# assert ("default", 1) in comp.predictors.keys()

83
miplearn/extractors.py
View File

@@ -1,6 +1,6 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
import numpy as np
from abc import ABC, abstractmethod
@@ -24,9 +24,26 @@ class Extractor(ABC):
result[category] = []
result[category] += [(var, index)]
return result
@staticmethod
def merge(partial_results, vertical=False):
results = {}
all_categories = set()
for pr in partial_results:
all_categories |= pr.keys()
for category in all_categories:
results[category] = []
for pr in partial_results:
if category in pr.keys():
results[category] += [pr[category]]
if vertical:
results[category] = np.vstack(results[category])
else:
results[category] = np.hstack(results[category])
return results
class UserFeaturesExtractor(Extractor):
class VariableFeaturesExtractor(Extractor):
def extract(self,
instances,
models=None,
@@ -45,6 +62,7 @@ class UserFeaturesExtractor(Extractor):
result[category] += [np.hstack([
instance_features,
instance.get_variable_features(var, index),
instance.lp_solution[str(var)][index],
])]
for category in result.keys():
result[category] = np.vstack(result[category])
@@ -52,8 +70,13 @@ class UserFeaturesExtractor(Extractor):
class SolutionExtractor(Extractor):
def extract(self, instances, models):
def __init__(self, relaxation=False):
self.relaxation = relaxation
def extract(self, instances, models=None):
result = {}
if models is None:
models = [instance.to_model() for instance in instances]
for (index, instance) in enumerate(instances):
model = models[index]
var_split = self.split_variables(instance, model)
@@ -61,10 +84,48 @@ class SolutionExtractor(Extractor):
if category not in result.keys():
result[category] = []
for (var, index) in var_index_pairs:
result[category] += [[
1 - var[index].value,
var[index].value,
]]
if self.relaxation:
v = instance.lp_solution[str(var)][index]
else:
v = instance.solution[str(var)][index]
if v is None:
result[category] += [[0, 0]]
else:
result[category] += [[1 - v, v]]
for category in result.keys():
result[category] = np.vstack(result[category])
return result
return result
class CombinedExtractor(Extractor):
def __init__(self, extractors):
self.extractors = extractors
def extract(self, instances, models):
return self.merge([ex.extract(instances, models)
for ex in self.extractors])
class InstanceFeaturesExtractor(Extractor):
def extract(self, instances, models=None):
return np.vstack([
np.hstack([
instance.get_instance_features(),
instance.lp_value,
])
for instance in instances
])
class ObjectiveValueExtractor(Extractor):
def __init__(self, kind="lp"):
assert kind in ["lower bound", "upper bound", "lp"]
self.kind = kind
def extract(self, instances, models=None):
if self.kind == "lower bound":
return np.array([[instance.lower_bound] for instance in instances])
if self.kind == "upper bound":
return np.array([[instance.upper_bound] for instance in instances])
if self.kind == "lp":
return np.array([[instance.lp_value] for instance in instances])

6
miplearn/instance.py
View File

@@ -1,6 +1,6 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from abc import ABC, abstractmethod

6
miplearn/problems/__init__.py
View File

@@ -1,3 +1,3 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

28
miplearn/problems/knapsack.py
View File

@@ -1,6 +1,6 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
import miplearn
from miplearn import Instance
@@ -13,7 +13,7 @@ from scipy.stats.distributions import rv_frozen
class ChallengeA:
"""
- 250 variables, 10 constraints, fixed weights
- w ~ U(100, 900), jitter ~ U(-100, 100)
- w ~ U(0, 1000), jitter ~ U(0.95, 1.05)
- K = 500, u ~ U(0., 1.)
- alpha = 0.25
"""
@@ -25,12 +25,12 @@ class ChallengeA:
np.random.seed(seed)
self.gen = MultiKnapsackGenerator(n=randint(low=250, high=251),
m=randint(low=10, high=11),
w=uniform(loc=100.0, scale=900.0),
w=uniform(loc=0.0, scale=1000.0),
K=uniform(loc=500.0, scale=0.0),
u=uniform(loc=0.0, scale=1.0),
alpha=uniform(loc=0.25, scale=0.0),
fix_w=True,
w_jitter=uniform(loc=-100.0, scale=200.0),
w_jitter=uniform(loc=0.95, scale=0.1),
)
np.random.seed(seed + 1)
self.training_instances = self.gen.generate(n_training_instances)
@@ -104,7 +104,7 @@ class MultiKnapsackGenerator:
u=uniform(loc=0.0, scale=1.0),
alpha=uniform(loc=0.25, scale=0.0),
fix_w=False,
w_jitter=randint(low=0, high=1),
w_jitter=uniform(loc=1.0, scale=0.0),
round=True,
):
"""Initialize the problem generator.
@@ -133,7 +133,7 @@ class MultiKnapsackGenerator:
be completely identical.
If a probability distribution w_jitter is provided, then item weights will be set to
w[i,j] + gamma[i,j] where gamma[i,j] is sampled from w_jitter. When combined with
w[i,j] * gamma[i,j] where gamma[i,j] is sampled from w_jitter. When combined with
fix_w=True, this argument may be used to generate instances where the weight of each item
is roughly the same, but not exactly identical, across all instances. The prices of the
items and the capacities of the knapsacks will be calculated as above, but using these
@@ -186,10 +186,14 @@ class MultiKnapsackGenerator:
self.fix_n = self.n.rvs()
self.fix_m = self.m.rvs()
self.fix_w = np.array([self.w.rvs(self.fix_n) for _ in range(self.fix_m)])
self.fix_u = self.u.rvs(self.fix_n)
self.fix_K = self.K.rvs()
else:
self.fix_n = None
self.fix_m = None
self.fix_w = None
self.fix_u = None
self.fix_K = None
def generate(self, n_samples):
def _sample():
@@ -197,13 +201,15 @@ class MultiKnapsackGenerator:
n = self.fix_n
m = self.fix_m
w = self.fix_w
u = self.fix_u
K = self.fix_K
else:
n = self.n.rvs()
m = self.m.rvs()
w = np.array([self.w.rvs(n) for _ in range(m)])
w = w + np.array([self.w_jitter.rvs(n) for _ in range(m)])
K = self.K.rvs()
u = self.u.rvs(n)
u = self.u.rvs(n)
K = self.K.rvs()
w = w * np.array([self.w_jitter.rvs(n) for _ in range(m)])
alpha = self.alpha.rvs(m)
p = np.array([w[:,j].sum() / m + K * u[j] for j in range(n)])
b = np.array([w[i,:].sum() * alpha[i] for i in range(m)])

6
miplearn/problems/stab.py
View File

@@ -1,6 +1,6 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
import numpy as np
import pyomo.environ as pe

4
miplearn/problems/tests/__init__.py
View File

@@ -0,0 +1,4 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

35
miplearn/problems/tests/test_knapsack.py
View File

@@ -1,6 +1,6 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import LearningSolver
from miplearn.problems.knapsack import MultiKnapsackGenerator, MultiKnapsackInstance
@@ -23,32 +23,3 @@ def test_knapsack_generator():
assert round(np.mean(w_sum), -1) == 500.
assert round(np.mean(p_sum), -1) == 1250.
assert round(np.mean(b_sum), -3) == 25000.
def test_knapsack_fixed_weights_jitter():
gen = MultiKnapsackGenerator(n=randint(low=50, high=51),
m=randint(low=10, high=11),
w=randint(low=0, high=1000),
K=randint(low=500, high=501),
u=uniform(loc=1.0, scale=0.0),
alpha=uniform(loc=0.50, scale=0.0),
fix_w=True,
w_jitter=randint(low=0, high=1),
)
instances = gen.generate(100)
w = [instance.weights[0,0] for instance in instances]
assert np.std(w) == 0.
gen = MultiKnapsackGenerator(n=randint(low=1, high=2),
m=randint(low=10, high=11),
w=randint(low=1000, high=1001),
K=randint(low=500, high=501),
u=uniform(loc=1.0, scale=0.0),
alpha=uniform(loc=0.50, scale=0.0),
fix_w=True,
w_jitter=randint(low=0, high=1001),
)
instances = gen.generate(5_000)
w = [instance.weights[0,0] for instance in instances]
assert round(np.std(w), -1) == 290.
assert round(np.mean(w), -2) == 1500.

6
miplearn/problems/tests/test_stab.py
View File

@@ -1,6 +1,6 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import LearningSolver
from miplearn.problems.stab import MaxWeightStableSetInstance

299
miplearn/solvers.py
View File

@@ -1,23 +1,170 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from .transformers import PerVariableTransformer
from .warmstart import WarmStartComponent
from .branching import BranchPriorityComponent
from . import ObjectiveValueComponent, PrimalSolutionComponent
import pyomo.environ as pe
import numpy as np
from pyomo.core import Var
from copy import deepcopy
import pickle
from scipy.stats import randint
from p_tqdm import p_map
import logging
logger = logging.getLogger(__name__)
def _gurobi_factory():
solver = pe.SolverFactory('gurobi_persistent')
solver.options["threads"] = 4
solver.options["Seed"] = randint(low=0, high=1000).rvs()
return solver
class InternalSolver:
def __init__(self):
self.is_warm_start_available = False
self.model = None
pass
def solve_lp(self, tee=False):
# Relax domain
from pyomo.core.base.set_types import Reals
original_domain = {}
for var in self.model.component_data_objects(Var):
original_domain[str(var)] = var.domain
lb, ub = var.bounds
var.setlb(lb)
var.setub(ub)
var.domain = Reals
# Solve LP relaxation
self.solver.set_instance(self.model)
results = self.solver.solve(tee=tee)
# Restore domains
for var in self.model.component_data_objects(Var):
var.domain = original_domain[str(var)]
# Reload original model
self.solver.set_instance(self.model)
return {
"Optimal value": results["Problem"][0]["Lower bound"],
}
def clear_values(self):
for var in self.model.component_objects(Var):
for index in var:
var[index].value = None
def get_solution(self):
solution = {}
for var in self.model.component_objects(Var):
solution[str(var)] = {}
for index in var:
solution[str(var)][index] = var[index].value
return solution
def set_warm_start(self, ws):
self.is_warm_start_available = True
self.clear_values()
count_total, count_fixed = 0, 0
for var in ws.keys():
for index in var:
count_total += 1
var[index].value = ws[var][index]
if ws[var][index] is not None:
count_fixed += 1
logger.info("Setting start values for %d variables (out of %d)" %
(count_fixed, count_total))
def set_model(self, model):
self.model = model
self.solver.set_instance(model)
def fix(self, ws):
count_total, count_fixed = 0, 0
for var in ws.keys():
for index in var:
count_total += 1
if ws[var][index] is None:
continue
count_fixed += 1
var[index].fix(ws[var][index])
self.solver.update_var(var[index])
logger.info("Fixing values for %d variables (out of %d)" %
(count_fixed, count_total))
class GurobiSolver(InternalSolver):
def __init__(self):
super().__init__()
self.solver = pe.SolverFactory('gurobi_persistent')
self.solver.options["Seed"] = randint(low=0, high=1000).rvs()
def set_threads(self, threads):
self.solver.options["Threads"] = threads
def set_time_limit(self, time_limit):
self.solver.options["TimeLimit"] = time_limit
def set_gap_tolerance(self, gap_tolerance):
self.solver.options["MIPGap"] = gap_tolerance
def solve(self, tee=False):
results = self.solver.solve(tee=tee, warmstart=self.is_warm_start_available)
return {
"Lower bound": results["Problem"][0]["Lower bound"],
"Upper bound": results["Problem"][0]["Upper bound"],
"Wallclock time": results["Solver"][0]["Wallclock time"],
"Nodes": self.solver._solver_model.getAttr("NodeCount"),
}
def _load_vars(self):
var_map = self._pyomo_var_to_solver_var_map
ref_vars = self._referenced_variables
vars_to_load = var_map.keys()
gurobi_vars_to_load = [var_map[pyomo_var] for pyomo_var in vars_to_load]
vals = self._solver_model.getAttr("X", gurobi_vars_to_load)
for var, val in zip(vars_to_load, vals):
if ref_vars[var] > 0:
var.stale = False
var.value = val
class CPLEXSolver(InternalSolver):
def __init__(self):
super().__init__()
import cplex
self.solver = pe.SolverFactory('cplex_persistent')
self.solver.options["randomseed"] = randint(low=0, high=1000).rvs()
def set_threads(self, threads):
self.solver.options["threads"] = threads
def set_time_limit(self, time_limit):
self.solver.options["timelimit"] = time_limit
def set_gap_tolerance(self, gap_tolerance):
self.solver.options["mip_tolerances_mipgap"] = gap_tolerance
def solve(self, tee=False):
results = self.solver.solve(tee=tee, warmstart=self.is_warm_start_available)
return {
"Lower bound": results["Problem"][0]["Lower bound"],
"Upper bound": results["Problem"][0]["Upper bound"],
"Wallclock time": results["Solver"][0]["Wallclock time"],
"Nodes": 1,
}
def solve_lp(self, tee=False):
import cplex
lp = self.solver._solver_model
var_types = lp.variables.get_types()
n_vars = len(var_types)
lp.set_problem_type(cplex.Cplex.problem_type.LP)
results = self.solver.solve(tee=tee)
lp.variables.set_types(zip(range(n_vars), var_types))
return {
"Optimal value": results["Problem"][0]["Lower bound"],
}
class LearningSolver:
"""
@@ -26,64 +173,92 @@ class LearningSolver:
"""
def __init__(self,
threads=None,
time_limit=None,
gap_limit=None,
internal_solver_factory=_gurobi_factory,
components=None,
mode=None):
gap_tolerance=None,
mode="exact",
solver="gurobi",
threads=4,
time_limit=None,
):
self.is_persistent = None
self.internal_solver = None
self.components = components
self.internal_solver_factory = internal_solver_factory
self.mode = mode
self.internal_solver = None
self.internal_solver_factory = solver
self.threads = threads
self.time_limit = time_limit
self.gap_limit = gap_limit
self.gap_tolerance = gap_tolerance
self.tee = False
self.training_instances = []
if self.components is not None:
assert isinstance(self.components, dict)
else:
self.components = {
"warm-start": WarmStartComponent(),
#"branch-priority": BranchPriorityComponent(),
"ObjectiveValue": ObjectiveValueComponent(),
"PrimalSolution": PrimalSolutionComponent(),
}
if mode is not None:
assert mode in ["exact", "heuristic"]
for component in self.components.values():
component.mode = mode
assert self.mode in ["exact", "heuristic"]
for component in self.components.values():
component.mode = self.mode
def _create_solver(self):
self.internal_solver = self.internal_solver_factory()
self.is_persistent = hasattr(self.internal_solver, "set_instance")
if self.threads is not None:
self.internal_solver.options["Threads"] = self.threads
def _create_internal_solver(self):
if self.internal_solver_factory == "cplex":
solver = CPLEXSolver()
elif self.internal_solver_factory == "gurobi":
solver = GurobiSolver()
else:
raise Exception("solver %s not supported" % solver_factory)
solver.set_threads(self.threads)
if self.time_limit is not None:
self.internal_solver.options["TimeLimit"] = self.time_limit
if self.gap_limit is not None:
self.internal_solver.options["MIPGap"] = self.gap_limit
solver.set_time_limit(self.time_limit)
if self.gap_tolerance is not None:
solver.set_gap_tolerance(self.gap_tolerance)
return solver
def solve(self, instance, tee=False):
model = instance.to_model()
self._create_solver()
if self.is_persistent:
self.internal_solver.set_instance(model)
def solve(self,
instance,
model=None,
tee=False,
relaxation_only=False,
):
if model is None:
model = instance.to_model()
self.tee = tee
self.internal_solver = self._create_internal_solver()
self.internal_solver.set_model(model)
# Solve LP relaxation
results = self.internal_solver.solve_lp(tee=tee)
instance.lp_solution = self.internal_solver.get_solution()
instance.lp_value = results["Optimal value"]
# Invoke before_solve callbacks
for component in self.components.values():
component.before_solve(self, instance, model)
if self.is_persistent:
solve_results = self.internal_solver.solve(tee=tee, warmstart=True)
else:
solve_results = self.internal_solver.solve(model, tee=tee, warmstart=True)
if relaxation_only:
return results
solve_results["Solver"][0]["Nodes"] = self.internal_solver._solver_model.getAttr("NodeCount")
# Solver original MIP
results = self.internal_solver.solve(tee=tee)
# Read MIP solution and bounds
instance.lower_bound = results["Lower bound"]
instance.upper_bound = results["Upper bound"]
instance.solution = self.internal_solver.get_solution()
# Invoke after_solve callbacks
for component in self.components.values():
component.after_solve(self, instance, model)
# Store instance for future training
self.training_instances += [instance]
return solve_results
return results
def parallel_solve(self,
instances,
@@ -99,11 +274,26 @@ class LearningSolver:
solver.internal_solver = None
if not collect_training_data:
solver.components = {}
return solver, results
return {
"Solver": solver,
"Results": results,
"Solution": instance.solution,
"LP solution": instance.lp_solution,
"LP value": instance.lp_value,
"Upper bound": instance.upper_bound,
"Lower bound": instance.lower_bound,
}
solver_result_pairs = p_map(_process, instances, num_cpus=n_jobs, desc=label)
subsolvers = [p[0] for p in solver_result_pairs]
results = [p[1] for p in solver_result_pairs]
p_map_results = p_map(_process, instances, num_cpus=n_jobs, desc=label)
subsolvers = [p["Solver"] for p in p_map_results]
results = [p["Results"] for p in p_map_results]
for (idx, r) in enumerate(p_map_results):
instances[idx].solution = r["Solution"]
instances[idx].lp_solution = r["LP solution"]
instances[idx].lp_value = r["LP value"]
instances[idx].lower_bound = r["Lower bound"]
instances[idx].upper_bound = r["Upper bound"]
for (name, component) in self.components.items():
subcomponents = [subsolver.components[name]
@@ -113,9 +303,13 @@ class LearningSolver:
return results
def fit(self, n_jobs=1):
def fit(self, training_instances=None):
if training_instances is None:
training_instances = self.training_instances
if len(training_instances) == 0:
return
for component in self.components.values():
component.fit(self, n_jobs=n_jobs)
component.fit(training_instances)
def save_state(self, filename):
with open(filename, "wb") as file:
@@ -133,3 +327,4 @@ class LearningSolver:
continue
else:
self.components[component_name].merge([component])

4
miplearn/tests/__init__.py
View File

@@ -0,0 +1,4 @@
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

11
miplearn/tests/test_benchmark.py
View File

@@ -1,9 +1,8 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import LearningSolver, BenchmarkRunner
from miplearn.warmstart import KnnWarmStartPredictor
from miplearn.problems.stab import MaxWeightStableSetGenerator
from scipy.stats import randint
import numpy as np
@@ -30,11 +29,11 @@ def test_benchmark():
benchmark = BenchmarkRunner(test_solvers)
benchmark.load_state("data.bin")
benchmark.parallel_solve(test_instances, n_jobs=2, n_trials=2)
assert benchmark.raw_results().values.shape == (12,12)
assert benchmark.raw_results().values.shape == (12,13)
benchmark.save_results("/tmp/benchmark.csv")
assert os.path.isfile("/tmp/benchmark.csv")
benchmark = BenchmarkRunner(test_solvers)
benchmark.load_results("/tmp/benchmark.csv")
assert benchmark.raw_results().values.shape == (12,12)
assert benchmark.raw_results().values.shape == (12,13)

49
miplearn/tests/test_branching.py
View File

@@ -1,49 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
from miplearn import BranchPriorityComponent, LearningSolver
from miplearn.problems.knapsack import KnapsackInstance
import numpy as np
import tempfile
def _get_instances():
return [
KnapsackInstance(
weights=[23., 26., 20., 18.],
prices=[505., 352., 458., 220.],
capacity=67.,
),
] * 2
def test_branching():
instances = _get_instances()
component = BranchPriorityComponent()
for instance in instances:
component.after_solve(None, instance, None)
component.fit(None)
for key in ["default"]:
assert key in component.x_train.keys()
assert key in component.y_train.keys()
assert component.x_train[key].shape == (8, 4)
assert component.y_train[key].shape == (8, 1)
def test_branch_priority_save_load():
state_file = tempfile.NamedTemporaryFile(mode="r")
solver = LearningSolver(components={"branch-priority": BranchPriorityComponent()})
solver.parallel_solve(_get_instances(), n_jobs=2)
solver.fit()
comp = solver.components["branch-priority"]
assert comp.x_train["default"].shape == (8, 4)
assert comp.y_train["default"].shape == (8, 1)
assert "default" in comp.predictors.keys()
solver.save_state(state_file.name)
solver = LearningSolver(components={"branch-priority": BranchPriorityComponent()})
solver.load_state(state_file.name)
comp = solver.components["branch-priority"]
assert comp.x_train["default"].shape == (8, 4)
assert comp.y_train["default"].shape == (8, 1)
assert "default" in comp.predictors.keys()

70
miplearn/tests/test_extractors.py
View File

@@ -1,16 +1,20 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn.problems.knapsack import KnapsackInstance
from miplearn import (UserFeaturesExtractor,
SolutionExtractor)
from miplearn import (LearningSolver,
SolutionExtractor,
CombinedExtractor,
InstanceFeaturesExtractor,
VariableFeaturesExtractor,
)
import numpy as np
import pyomo.environ as pe
def _get_instances():
return [
instances = [
KnapsackInstance(weights=[1., 2., 3.],
prices=[10., 20., 30.],
capacity=2.5,
@@ -20,26 +24,16 @@ def _get_instances():
capacity=4.5,
),
]
def test_user_features():
instances = _get_instances()
extractor = UserFeaturesExtractor()
features = extractor.extract(instances)
assert isinstance(features, dict)
assert "default" in features.keys()
assert isinstance(features["default"], np.ndarray)
assert features["default"].shape == (6, 4)
def test_solution_extractor():
instances = _get_instances()
models = [instance.to_model() for instance in instances]
for model in models:
solver = pe.SolverFactory("cbc")
solver.solve(model)
extractor = SolutionExtractor()
features = extractor.extract(instances, models)
solver = LearningSolver()
for (i, instance) in enumerate(instances):
solver.solve(instances[i], models[i])
return instances, models
def test_solution_extractor():
instances, models = _get_instances()
features = SolutionExtractor().extract(instances, models)
assert isinstance(features, dict)
assert "default" in features.keys()
assert isinstance(features["default"], np.ndarray)
@@ -52,3 +46,29 @@ def test_solution_extractor():
0., 1.,
1., 0.,
]
def test_combined_extractor():
instances, models = _get_instances()
extractor = CombinedExtractor(extractors=[VariableFeaturesExtractor(),
SolutionExtractor()])
features = extractor.extract(instances, models)
assert isinstance(features, dict)
assert "default" in features.keys()
assert isinstance(features["default"], np.ndarray)
assert features["default"].shape == (6, 7)
def test_instance_features_extractor():
instances, models = _get_instances()
features = InstanceFeaturesExtractor().extract(instances)
assert features.shape == (2,3)
def test_variable_features_extractor():
instances, models = _get_instances()
features = VariableFeaturesExtractor().extract(instances)
assert isinstance(features, dict)
assert "default" in features
assert features["default"].shape == (6,5)

93
miplearn/tests/test_solver.py
View File

@@ -1,12 +1,9 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
# MIPLearn: Extensible Framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.
from miplearn import LearningSolver
from miplearn import LearningSolver, BranchPriorityComponent
from miplearn.problems.knapsack import KnapsackInstance
from miplearn.branching import BranchPriorityComponent
from miplearn.warmstart import WarmStartComponent
import numpy as np
def _get_instance():
@@ -16,48 +13,62 @@ def _get_instance():
capacity=67.,
)
def test_solver():
instance = _get_instance()
solver = LearningSolver()
solver.solve(instance)
solver.fit()
solver.solve(instance)
for mode in ["exact", "heuristic"]:
for internal_solver in ["cplex", "gurobi"]:
solver = LearningSolver(time_limit=300,
gap_tolerance=1e-3,
threads=1,
solver=internal_solver,
mode=mode,
)
results = solver.solve(instance)
assert instance.solution["x"][0] == 1.0
assert instance.solution["x"][1] == 0.0
assert instance.solution["x"][2] == 1.0
assert instance.solution["x"][3] == 1.0
assert instance.lower_bound == 1183.0
assert instance.upper_bound == 1183.0
def test_solve_save_load_state():
instance = _get_instance()
components_before = {
"warm-start": WarmStartComponent(),
"branch-priority": BranchPriorityComponent(),
}
solver = LearningSolver(components=components_before)
solver.solve(instance)
solver.fit()
solver.save_state("/tmp/knapsack_train.bin")
prev_x_train_len = len(solver.components["warm-start"].x_train)
prev_y_train_len = len(solver.components["warm-start"].y_train)
assert round(instance.lp_solution["x"][0], 3) == 1.000
assert round(instance.lp_solution["x"][1], 3) == 0.923
assert round(instance.lp_solution["x"][2], 3) == 1.000
assert round(instance.lp_solution["x"][3], 3) == 0.000
assert round(instance.lp_value, 3) == 1287.923
solver.fit()
solver.solve(instance)
# def test_solve_save_load_state():
# instance = _get_instance()
# components_before = {
# "warm-start": WarmStartComponent(),
# }
# solver = LearningSolver(components=components_before)
# solver.solve(instance)
# solver.fit()
# solver.save_state("/tmp/knapsack_train.bin")
# prev_x_train_len = len(solver.components["warm-start"].x_train)
# prev_y_train_len = len(solver.components["warm-start"].y_train)
components_after = {
"warm-start": WarmStartComponent(),
}
solver = LearningSolver(components=components_after)
solver.load_state("/tmp/knapsack_train.bin")
assert len(solver.components.keys()) == 1
assert len(solver.components["warm-start"].x_train) == prev_x_train_len
assert len(solver.components["warm-start"].y_train) == prev_y_train_len
# components_after = {
# "warm-start": WarmStartComponent(),
# }
# solver = LearningSolver(components=components_after)
# solver.load_state("/tmp/knapsack_train.bin")
# assert len(solver.components.keys()) == 1
# assert len(solver.components["warm-start"].x_train) == prev_x_train_len
# assert len(solver.components["warm-start"].y_train) == prev_y_train_len
def test_parallel_solve():
instances = [_get_instance() for _ in range(10)]
solver = LearningSolver()
results = solver.parallel_solve(instances, n_jobs=3)
assert len(results) == 10
assert len(solver.components["warm-start"].x_train["default"]) == 40
assert len(solver.components["warm-start"].y_train["default"]) == 40
for instance in instances:
assert len(instance.solution["x"].keys()) == 4
def test_solver_random_branch_priority():
instance = _get_instance()
components = {
"warm-start": BranchPriorityComponent(),
}
solver = LearningSolver(components=components)
solver.solve(instance)
solver.fit()

37
miplearn/tests/test_warmstart.py
View File

@@ -1,37 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
from miplearn import WarmStartComponent, LearningSolver
from miplearn.problems.knapsack import KnapsackInstance
import numpy as np
import tempfile
def _get_instances():
return [
KnapsackInstance(
weights=[23., 26., 20., 18.],
prices=[505., 352., 458., 220.],
capacity=67.,
),
] * 2
def test_warm_start_save_load():
state_file = tempfile.NamedTemporaryFile(mode="r")
solver = LearningSolver(components={"warm-start": WarmStartComponent()})
solver.parallel_solve(_get_instances(), n_jobs=2)
solver.fit()
comp = solver.components["warm-start"]
assert comp.x_train["default"].shape == (8, 4)
assert comp.y_train["default"].shape == (8, 2)
assert "default" in comp.predictors.keys()
solver.save_state(state_file.name)
solver = LearningSolver(components={"warm-start": WarmStartComponent()})
solver.load_state(state_file.name)
comp = solver.components["warm-start"]
assert comp.x_train["default"].shape == (8, 4)
assert comp.y_train["default"].shape == (8, 2)
assert "default" in comp.predictors.keys()

88
miplearn/tests/test_warmstart_knn.py
View File

@@ -1,88 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
from miplearn.warmstart import KnnWarmStartPredictor
from sklearn.metrics import accuracy_score, precision_score
import numpy as np
def test_knn_with_consensus():
x_train = np.array([
[0.0, 0.0],
[0.1, 0.0],
[0.0, 0.1],
[1.0, 1.0],
])
y_train = np.array([
[0., 1.],
[0., 1.],
[0., 1.],
[1., 0.],
])
ws = KnnWarmStartPredictor(k=3, thr_clip=[0.75, 0.75])
ws.fit(x_train, y_train)
x_test = np.array([[0.0, 0.0]])
y_test = np.array([[0, 1]])
assert (ws.predict(x_test) == y_test).all()
def test_knn_without_consensus():
x_train = np.array([
[0.0, 0.0],
[0.1, 0.1],
[0.9, 0.9],
[1.0, 1.0],
])
y_train = np.array([
[0., 1.],
[0., 1.],
[1., 0.],
[1., 0.],
])
ws = KnnWarmStartPredictor(k=4, thr_clip=[0.75, 0.75])
ws.fit(x_train, y_train)
x_test = np.array([[0.5, 0.5]])
y_test = np.array([[0, 0]])
assert (ws.predict(x_test) == y_test).all()
def test_knn_always_true():
x_train = np.array([
[0.0, 0.0],
[0.1, 0.1],
[0.9, 0.9],
[1.0, 1.0],
])
y_train = np.array([
[1., 0.],
[1., 0.],
[1., 0.],
[1., 0.],
])
ws = KnnWarmStartPredictor(k=4, thr_clip=[0.75, 0.75])
ws.fit(x_train, y_train)
x_test = np.array([[0.5, 0.5]])
y_test = np.array([[1, 0]])
assert (ws.predict(x_test) == y_test).all()
def test_knn_always_false():
x_train = np.array([
[0.0, 0.0],
[0.1, 0.1],
[0.9, 0.9],
[1.0, 1.0],
])
y_train = np.array([
[0., 1.],
[0., 1.],
[0., 1.],
[0., 1.],
])
ws = KnnWarmStartPredictor(k=4, thr_clip=[0.75, 0.75])
ws.fit(x_train, y_train)
x_test = np.array([[0.5, 0.5]])
y_test = np.array([[0, 1]])
assert (ws.predict(x_test) == y_test).all()

64
miplearn/tests/test_warmstart_logistic.py
View File

@@ -1,64 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
from miplearn.warmstart import LogisticWarmStartPredictor
from sklearn.metrics import accuracy_score, precision_score
import numpy as np
def _generate_dataset(ground_truth, n_samples=10_000):
x_train = np.random.rand(n_samples,5)
x_test = np.random.rand(n_samples,5)
y_train = ground_truth(x_train)
y_test = ground_truth(x_test)
return x_train, y_train, x_test, y_test
def _is_sum_greater_than_two(x):
y = (np.sum(x, axis=1) > 2.0).astype(int)
return np.vstack([y, 1 - y]).transpose()
def _always_zero(x):
y = np.zeros((1, x.shape[0]))
return np.vstack([y, 1 - y]).transpose()
def _random_values(x):
y = np.random.randint(2, size=x.shape[0])
return np.vstack([y, 1 - y]).transpose()
def test_logistic_ws_with_balanced_labels():
x_train, y_train, x_test, y_test = _generate_dataset(_is_sum_greater_than_two)
ws = LogisticWarmStartPredictor()
ws.fit(x_train, y_train)
y_pred = ws.predict(x_test)
assert accuracy_score(y_test[:,0], y_pred[:,0]) > 0.99
assert accuracy_score(y_test[:,1], y_pred[:,1]) > 0.99
def test_logistic_ws_with_unbalanced_labels():
x_train, y_train, x_test, y_test = _generate_dataset(_always_zero)
ws = LogisticWarmStartPredictor()
ws.fit(x_train, y_train)
y_pred = ws.predict(x_test)
assert accuracy_score(y_test[:,0], y_pred[:,0]) == 1.0
assert accuracy_score(y_test[:,1], y_pred[:,1]) == 1.0
def test_logistic_ws_with_unpredictable_labels():
x_train, y_train, x_test, y_test = _generate_dataset(_random_values)
ws = LogisticWarmStartPredictor()
ws.fit(x_train, y_train)
y_pred = ws.predict(x_test)
assert np.sum(y_pred) == 0
def test_logistic_ws_with_small_sample_size():
x_train, y_train, x_test, y_test = _generate_dataset(_random_values, n_samples=3)
ws = LogisticWarmStartPredictor()
ws.fit(x_train, y_train)
y_pred = ws.predict(x_test)
assert np.sum(y_pred) == 0

65
miplearn/transformers.py
View File

@@ -1,65 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
import numpy as np
from pyomo.core import Var
class PerVariableTransformer:
"""
Class that converts a miplearn.Instance into a matrix of features that is suitable
for training machine learning models that make one decision per decision variable.
"""
def __init__(self):
pass
def transform_instance(self, instance, var_index_pairs):
instance_features = self._get_instance_features(instance)
variable_features = self._get_variable_features(instance, var_index_pairs)
return np.vstack([
np.hstack([instance_features, vf])
for vf in variable_features
])
@staticmethod
def _get_instance_features(instance):
features = instance.get_instance_features()
assert isinstance(features, np.ndarray)
return features
@staticmethod
def _get_variable_features(instance, var_index_pairs):
features = []
expected_shape = None
for (var, index) in var_index_pairs:
vf = instance.get_variable_features(var, index)
assert isinstance(vf, np.ndarray)
if expected_shape is None:
assert len(vf.shape) == 1
expected_shape = vf.shape
else:
assert vf.shape == expected_shape
features += [vf]
return np.array(features)
@staticmethod
def transform_solution(var_index_pairs):
y = []
for (var, index) in var_index_pairs:
y += [[1 - var[index].value, var[index].value]]
return np.array(y)
@staticmethod
def split_variables(instance, model):
result = {}
for var in model.component_objects(Var):
for index in var:
category = instance.get_variable_category(var, index)
if category is None:
continue
if category not in result.keys():
result[category] = []
result[category] += [(var, index)]
return result

211
miplearn/warmstart.py
View File

@@ -1,211 +0,0 @@
# MIPLearn, an extensible framework for Learning-Enhanced Mixed-Integer Optimization
# Copyright (C) 2019-2020 Argonne National Laboratory. All rights reserved.
# Written by Alinson S. Xavier <axavier@anl.gov>
from . import Component
from .transformers import PerVariableTransformer
from abc import ABC, abstractmethod
from copy import deepcopy
import numpy as np
from sklearn.pipeline import make_pipeline
from sklearn.linear_model import LogisticRegression
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import cross_val_score
from sklearn.neighbors import KNeighborsClassifier
from tqdm.auto import tqdm
class WarmStartPredictor(ABC):
def __init__(self, thr_clip=[0.50, 0.50]):
self.models = [None, None]
self.thr_clip = thr_clip
def fit(self, x_train, y_train):
assert isinstance(x_train, np.ndarray)
assert isinstance(y_train, np.ndarray)
y_train = y_train.astype(int)
assert y_train.shape[0] == x_train.shape[0]
assert y_train.shape[1] == 2
for i in [0,1]:
self.models[i] = self._fit(x_train, y_train[:, i], i)
def predict(self, x_test):
assert isinstance(x_test, np.ndarray)
y_pred = np.zeros((x_test.shape[0], 2))
for i in [0,1]:
if isinstance(self.models[i], int):
y_pred[:, i] = self.models[i]
else:
y = self.models[i].predict_proba(x_test)[:,1]
y[y < self.thr_clip[i]] = 0.
y[y > 0.] = 1.
y_pred[:, i] = y
return y_pred.astype(int)
@abstractmethod
def _fit(self, x_train, y_train, label):
pass
class LogisticWarmStartPredictor(WarmStartPredictor):
def __init__(self,
min_samples=100,
thr_fix=[0.99, 0.99],
thr_balance=[0.80, 0.80],
thr_alpha=[0.50, 0.50],
):
super().__init__()
self.min_samples = min_samples
self.thr_fix = thr_fix
self.thr_balance = thr_balance
self.thr_alpha = thr_alpha
def _fit(self, x_train, y_train, label):
y_train_avg = np.average(y_train)
# If number of samples is too small, don't predict anything.
if x_train.shape[0] < self.min_samples:
return 0
# If vast majority of observations are true, always return true.
if y_train_avg > self.thr_fix[label]:
return 1
# If dataset is not balanced enough, don't predict anything.
if y_train_avg < (1 - self.thr_balance[label]) or y_train_avg > self.thr_balance[label]:
return 0
reg = make_pipeline(StandardScaler(), LogisticRegression())
reg_score = np.mean(cross_val_score(reg, x_train, y_train, cv=5))
dummy_score = max(y_train_avg, 1 - y_train_avg)
reg_thr = 1. * self.thr_alpha[label] + dummy_score * (1 - self.thr_alpha[label])
# If cross-validation score is too low, don't predict anything.
if reg_score < reg_thr:
return 0
reg.fit(x_train, y_train.astype(int))
return reg
class KnnWarmStartPredictor(WarmStartPredictor):
def __init__(self, k=50,
thr_clip=[0.90, 0.90],
thr_fix=[0.99, 0.99],
):
super().__init__(thr_clip=thr_clip)
self.k = k
self.thr_fix = thr_fix
def _fit(self, x_train, y_train, label):
y_train_avg = np.average(y_train)
# If number of training samples is too small, don't predict anything.
if x_train.shape[0] < self.k:
return 0
# If vast majority of observations are true, always return true.
if y_train_avg > self.thr_fix[label]:
return 1
# If vast majority of observations are false, always return false.
if y_train_avg < (1 - self.thr_fix[label]):
return 0
knn = KNeighborsClassifier(n_neighbors=self.k)
knn.fit(x_train, y_train)
return knn
class WarmStartComponent(Component):
def __init__(self,
predictor_prototype=KnnWarmStartPredictor(),
mode="exact",
):
self.mode = mode
self.transformer = PerVariableTransformer()
self.x_train = {}
self.y_train = {}
self.predictors = {}
self.predictor_prototype = predictor_prototype
def before_solve(self, solver, instance, model):
var_split = self.transformer.split_variables(instance, model)
x_test = {}
# Collect training data (x_train) and build x_test
for category in var_split.keys():
var_index_pairs = var_split[category]
x = self.transformer.transform_instance(instance, var_index_pairs)
x_test[category] = x
if category not in self.x_train.keys():
self.x_train[category] = x
else:
assert x.shape[1] == self.x_train[category].shape[1]
self.x_train[category] = np.vstack([self.x_train[category], x])
# Predict solutions
for category in var_split.keys():
var_index_pairs = var_split[category]
if category in self.predictors.keys():
ws = self.predictors[category].predict(x_test[category])
assert ws.shape == (len(var_index_pairs), 2)
for i in range(len(var_index_pairs)):
var, index = var_index_pairs[i]
if self.mode == "heuristic":
if ws[i,0] == 1:
var[index].fix(0)
if solver.is_persistent:
solver.internal_solver.update_var(var[index])
elif ws[i,1] == 1:
var[index].fix(1)
if solver.is_persistent:
solver.internal_solver.update_var(var[index])
else:
if ws[i,0] == 1:
var[index].value = 0
elif ws[i,1] == 1:
var[index].value = 1
def after_solve(self, solver, instance, model):
var_split = self.transformer.split_variables(instance, model)
for category in var_split.keys():
var_index_pairs = var_split[category]
y = self.transformer.transform_solution(var_index_pairs)
if category not in self.y_train.keys():
self.y_train[category] = y
else:
self.y_train[category] = np.vstack([self.y_train[category], y])
def fit(self, solver, n_jobs=1):
for category in tqdm(self.x_train.keys(), desc="Warm start"):
x_train = self.x_train[category]
y_train = self.y_train[category]
self.predictors[category] = deepcopy(self.predictor_prototype)
self.predictors[category].fit(x_train, y_train)
def merge(self, other_components):
# Merge x_train and y_train
keys = set(self.x_train.keys())
for comp in other_components:
keys = keys.union(set(comp.x_train.keys()))
for key in keys:
x_train_submatrices = [comp.x_train[key]
for comp in other_components
if key in comp.x_train.keys()]
y_train_submatrices = [comp.y_train[key]
for comp in other_components
if key in comp.y_train.keys()]
if key in self.x_train.keys():
x_train_submatrices += [self.x_train[key]]
y_train_submatrices += [self.y_train[key]]
self.x_train[key] = np.vstack(x_train_submatrices)
self.y_train[key] = np.vstack(y_train_submatrices)
# Merge trained predictors
for comp in other_components:
for key in comp.predictors.keys():
if key not in self.predictors.keys():
self.predictors[key] = comp.predictors[key]