mlp with bounds to directly optimizer

constrained_linear_regression/constrained_multi_layer_perceptron.py

@@ -33,7 +33,6 @@ def fit(self, X, y, min_coef=None, max_coef=None):
         feature_count = X.shape[-1]
         self.min_coef_ = self._verify_coef(feature_count, min_coef, -np.inf).flatten()
         self.max_coef_ = self._verify_coef(feature_count, max_coef, np.inf).flatten()
-
         return self._constrained_fit(X, y, incremental=False)
 
     def _constrained_fit(self, X, y, incremental=False):
@@ -261,15 +260,35 @@ def _fit_constrained_stochastic(
 
         self._after_training()
 
-    def _update_coef_using_constrain(self, coef_grads):
-        for coef_idx, (min_coef, max_coef) in enumerate(
-            zip(self.min_coef_, self.max_coef_)
-        ):
-            # clipping is applied only to the first node.
-            coef_grads[0][coef_idx] = np.clip(
-                coef_grads[0][coef_idx], min_coef, max_coef
+    def _update_coef_using_constrain(self, packed_coef_inter=None):
+        if self.solver == "lbfgs":
+            begin = self._coef_indptr[0][0]
+            end = self._coef_indptr[0][1]
+            shape = self._coef_indptr[0][2]
+            unpacked_coef_inter = np.reshape(
+                packed_coef_inter[begin:end],
+                shape,
             )
 
+            for coef_idx, (min_coef, max_coef) in enumerate(
+                zip(self.min_coef_, self.max_coef_)
+            ):
+                # clipping is applied only to the first node.
+                unpacked_coef_inter[coef_idx] = np.clip(
+                    unpacked_coef_inter[coef_idx], min_coef, max_coef
+                )
+            packed_coef_inter[begin:end] = unpacked_coef_inter.reshape(-1)
+            return packed_coef_inter
+
+        else:
+            for coef_idx, (min_coef, max_coef) in enumerate(
+                zip(self.min_coef_, self.max_coef_)
+            ):
+                # clipping is applied only to the first node.
+                self.coefs_[0][coef_idx] = np.clip(
+                    self.coefs_[0][coef_idx], min_coef, max_coef
+                )
+
     def _after_training(self):
         # This is useful only for multi models.
         pass
@@ -304,11 +323,21 @@ def _fit_constrained_lbfgs(
         else:
             iprint = -1
 
+        bounds = [
+            [(min_val, max_val)] * self.hidden_layer_sizes[0]
+            for min_val, max_val in zip(self.min_coef_, self.max_coef_)
+        ]
+        bounds = list(chain(*bounds))
+        effective_bounds_len = len(bounds)
+        for _ in range(len(packed_coef_inter) - effective_bounds_len):
+            bounds.append((-np.inf, np.inf))
+
         opt_res = scipy.optimize.minimize(
             self._loss_grad_constrained_lbfgs,
             packed_coef_inter,
             method="L-BFGS-B",
             jac=True,
+            bounds=bounds,
             options={
                 "maxfun": self.max_fun,
                 "maxiter": self.max_iter,
@@ -369,7 +398,6 @@ def _loss_grad_constrained_lbfgs(
         loss, coef_grads, intercept_grads = self._backprop(
             X, y, activations, deltas, coef_grads, intercept_grads
         )
-        self._update_coef_using_constrain(coef_grads)  # it throws error.
         grad = _pack(coef_grads, intercept_grads)
         self._save_mse(loss)
         return loss, grad

constrained_linear_regression/multi_constrained_multi_layer_perceptron.py

@@ -1,8 +1,13 @@
 import numpy as np
+from itertools import chain
 
 from constrained_linear_regression.constrained_multi_layer_perceptron import (
     ConstrainedMultilayerPerceptron,
+    _pack,
 )
+import scipy.optimize
+
+from sklearn.utils.optimize import _check_optimize_result
 
 
 class MultiConstrainedMultilayerPerceptron(ConstrainedMultilayerPerceptron):
@@ -56,6 +61,71 @@ def _update_coef_using_constrain(self, coef_grads):
                 max_coef,
             )
 
+    def _fit_constrained_lbfgs(
+        self, X, y, activations, deltas, coef_grads, intercept_grads, layer_units
+    ):
+        # Store meta information for the parameters
+        self._coef_indptr = []
+        self._intercept_indptr = []
+        start = 0
+
+        # Save sizes and indices of coefficients for faster unpacking
+        for i in range(self.n_layers_ - 1):
+            n_fan_in, n_fan_out = layer_units[i], layer_units[i + 1]
+
+            end = start + (n_fan_in * n_fan_out)
+            self._coef_indptr.append((start, end, (n_fan_in, n_fan_out)))
+            start = end
+
+        # Save sizes and indices of intercepts for faster unpacking
+        for i in range(self.n_layers_ - 1):
+            end = start + layer_units[i + 1]
+            self._intercept_indptr.append((start, end))
+            start = end
+
+        # Run LBFGS
+        packed_coef_inter = _pack(self.coefs_, self.intercepts_)
+
+        if self.verbose is True or self.verbose >= 1:
+            iprint = 1
+        else:
+            iprint = -1
+
+        bounds = [
+            [(min_val, max_val)] * self.hidden_layer_sizes[0]
+            for min_val, max_val in zip(
+                self.min_coef_[
+                    MultiConstrainedMultilayerPerceptron.global_horizon_count
+                ],
+                self.max_coef_[
+                    MultiConstrainedMultilayerPerceptron.global_horizon_count
+                ],
+            )
+        ]
+        bounds = list(chain(*bounds))
+        effective_bounds_len = len(bounds)
+        for _ in range(len(packed_coef_inter) - effective_bounds_len):
+            bounds.append((-np.inf, np.inf))
+
+        opt_res = scipy.optimize.minimize(
+            self._loss_grad_constrained_lbfgs,
+            packed_coef_inter,
+            method="L-BFGS-B",
+            jac=True,
+            bounds=bounds,
+            options={
+                "maxfun": self.max_fun,
+                "maxiter": self.max_iter,
+                "iprint": iprint,
+                "gtol": self.tol,
+            },
+            args=(X, y, activations, deltas, coef_grads, intercept_grads),
+        )
+        self.n_iter_ = _check_optimize_result("lbfgs", opt_res, self.max_iter)
+        self.loss_ = opt_res.fun
+        self._unpack(opt_res.x)
+        self._after_training()
+
     def reset(self):
         MultiConstrainedMultilayerPerceptron.global_horizon_count = 0
         return f"horizon_count: {MultiConstrainedMultilayerPerceptron.global_horizon_count}"

tests/test_clr_basic.py

@@ -15,6 +15,9 @@
 from constrained_linear_regression.selective_drop_linear_regression import (
     SelectiveDropLinearRegression,
 )
+from constrained_linear_regression.selective_drop_multi_layer_perceptron import (
+    SelectiveDropMultilayerPerceptron,
+)
 
 
 def test_unconstrained():
@@ -29,7 +32,7 @@ def test_unconstrained():
     assert np.isclose(baseline.intercept_, model.intercept_)
 
 
-def test_nodrop():
+def test_nodrop_lr():
     X, Y = load_linnerud(return_X_y=True)
     y = Y[:, 0]
     model = SelectiveDropLinearRegression(nonnegative=False)
@@ -41,6 +44,56 @@ def test_nodrop():
     assert np.isclose(baseline.intercept_, model.intercept_)
 
 
+def test_nodrop_mlp_lbfgs(solver="lbfgs"):
+    X, Y = load_linnerud(return_X_y=True)
+    y = Y[:, 0]
+    random_state = 7
+    hidden_layer_sizes = (3,)
+    baseline = MLPRegressor(
+        solver=solver,
+        shuffle=False,
+        hidden_layer_sizes=hidden_layer_sizes,
+        random_state=random_state,
+    )
+    model = SelectiveDropMultilayerPerceptron(
+        solver=solver, hidden_layer_sizes=hidden_layer_sizes, random_state=random_state
+    )
+    model.fit(X, y)
+    baseline.fit(X, y)
+    for baseline_coef, model_coef in zip(baseline.coefs_, model.coefs_):
+        assert np.allclose(baseline_coef, model_coef)
+
+    for baseline_intercept, model_intercept in zip(
+        baseline.intercepts_, model.intercepts_
+    ):
+        assert np.allclose(baseline_intercept, model_intercept)
+
+
+def test_nodrop_mlp_sgd(solver="sgd"):
+    X, Y = load_linnerud(return_X_y=True)
+    y = Y[:, 0]
+    random_state = 7
+    hidden_layer_sizes = (3,)
+    baseline = MLPRegressor(
+        solver=solver,
+        shuffle=False,
+        hidden_layer_sizes=hidden_layer_sizes,
+        random_state=random_state,
+    )
+    model = SelectiveDropMultilayerPerceptron(
+        solver=solver, hidden_layer_sizes=hidden_layer_sizes, random_state=random_state
+    )
+    model.fit(X, y)
+    baseline.fit(X, y)
+    for baseline_coef, model_coef in zip(baseline.coefs_, model.coefs_):
+        assert np.allclose(baseline_coef, model_coef)
+
+    for baseline_intercept, model_intercept in zip(
+        baseline.intercepts_, model.intercepts_
+    ):
+        assert np.allclose(baseline_intercept, model_intercept)
+
+
 def test_positive():
     X, Y = load_linnerud(return_X_y=True)
     y = Y[:, 0]
@@ -112,6 +165,22 @@ def test_unconstrainedmlp_lbfgs(solver="lbfgs"):
         assert np.allclose(baseline_intercept, model_intercept)
 
 
+def test_constrainedmlp_lbfgs(solver="lbfgs"):
+    X, Y = load_linnerud(return_X_y=True)
+    y = Y[:, 0]
+    random_state = 7
+    hidden_layer_sizes = (2,)
+    min_coef = np.repeat(0, 3)
+    max_coef = np.repeat(2, 3)
+
+    max_coef[1] = 1
+    max_coef[2] = 0
+    model = ConstrainedMultilayerPerceptron(
+        solver=solver, hidden_layer_sizes=hidden_layer_sizes, random_state=random_state
+    )
+    model.fit(X, y, max_coef=max_coef, min_coef=min_coef)
+    assert np.all(model.coefs_[0] >= 0)
+
+
 if __name__ == "__main__":
-    test_unconstrainedmlp_lbfgs()
-    # test_unconstrainedmlp_sgd()
+    test_constrainedmlp_lbfgs()