From 33e5d1fefebc4c389bd97d7d3d88000bbd949824 Mon Sep 17 00:00:00 2001
From: Silvio Borer <silvio.borer@rohde-schwarz.com>
Date: Fri, 11 Sep 2020 11:12:25 +0200
Subject: [PATCH] add penalties/constraints for positive and negative functions

---
 pygam/penalties.py            | 54 +++++++++++++++++++++++++++++++++++
 pygam/tests/test_penalties.py | 31 ++++++++++++++++++++
 2 files changed, 85 insertions(+)

diff --git a/pygam/penalties.py b/pygam/penalties.py
index aee1d2a0..5ec4571c 100644
--- a/pygam/penalties.py
+++ b/pygam/penalties.py
@@ -211,6 +211,58 @@ def concave(n, coef):
     """
     return convexity_(n, coef, convex=False)
 
+def positive(n, coef):
+    """
+    Builds a penalty matrix for P-Splines with continuous features.
+    Penalizes violation of a positive feature function.
+
+    Parameters
+    ----------
+    n : int
+        number of splines
+    coef : array-like
+        coefficients of the feature function
+
+    Returns
+    -------
+    penalty matrix : sparse csc matrix of shape (n,n)
+    """
+    if n != len(coef.ravel()):
+        raise ValueError('dimension mismatch: expected n equals len(coef), '\
+                         'but found n = {}, coef.shape = {}.'\
+                         .format(n, coef.shape))
+    # only penalize the case where coef_i-1 < 0
+    mask = sp.sparse.diags((coef.ravel() < 0).astype(float))
+
+    D = sp.sparse.identity(n).tocsc() * mask
+    return D.dot(D.T).tocsc()
+
+def negative(n, coef):
+    """
+    Builds a penalty matrix for P-Splines with continuous features.
+    Penalizes violation of a negative feature function.
+
+    Parameters
+    ----------
+    n : int
+        number of splines
+    coef : array-like
+        coefficients of the feature function
+
+    Returns
+    -------
+    penalty matrix : sparse csc matrix of shape (n,n)
+    """
+    if n != len(coef.ravel()):
+        raise ValueError('dimension mismatch: expected n equals len(coef), '\
+                         'but found n = {}, coef.shape = {}.'\
+                         .format(n, coef.shape))
+    # only penalize the case where coef_i-1 > 0
+    mask = sp.sparse.diags((coef.ravel() > 0).astype(float))
+
+    D = sp.sparse.identity(n).tocsc() * mask
+    return D.dot(D.T).tocsc()
+
 # def circular(n, coef):
 #     """
 #     Builds a penalty matrix for P-Splines with continuous features.
@@ -341,5 +393,7 @@ def sparse_diff(array, n=1, axis=-1):
                'concave': concave,
                'monotonic_inc': monotonic_inc,
                'monotonic_dec': monotonic_dec,
+               'positive': positive,
+               'negative': negative,
                'none': none
               }
diff --git a/pygam/tests/test_penalties.py b/pygam/tests/test_penalties.py
index 68037619..8184911c 100644
--- a/pygam/tests/test_penalties.py
+++ b/pygam/tests/test_penalties.py
@@ -11,6 +11,8 @@
 from pygam.penalties import monotonic_dec
 from pygam.penalties import convex
 from pygam.penalties import concave
+from pygam.penalties import positive
+from pygam.penalties import negative
 from pygam.penalties import none
 from pygam.penalties import wrap_penalty
 
@@ -107,6 +109,35 @@ def test_concave(hepatitis_X_y):
     diffs = np.diff(Y, n=2)
     assert(((diffs <= 0) + np.isclose(diffs, 0.)).all())
 
+def test_positive(hepatitis_X_y):
+    """
+    check that positive constraint produces positive function
+    """
+    X, y = hepatitis_X_y
+
+    gam = LinearGAM(terms=s(0, constraints='positive'))
+    gam.fit(X, y)
+
+    XX = gam.generate_X_grid(term=0)
+    Y = gam.predict(np.sort(XX))
+    assert(((Y >= 0) + np.isclose(Y, 0.)).all())
+
+def test_negative(hepatitis_X_y):
+    """
+    check that negative constraint produces positive function
+    """
+    X, y = hepatitis_X_y
+    # invert sign to be able to fit negative function
+    X = -X
+    y = -y
+
+    gam = LinearGAM(terms=s(0, constraints='negative'))
+    gam.fit(X, y)
+
+    XX = gam.generate_X_grid(term=0)
+    Y = gam.predict(np.sort(XX))
+    assert(((Y <= 0) + np.isclose(Y, 0.)).all())
+
 
 # TODO penalties gives expected matrix structure
 # TODO circular constraints